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Abstract: We present new determinations of the MS charm quark mass using relativistic 
QCD sum rules at 0(ag) from the moments of the vector and the pseudoscalar current 
correlators. We use available experimental measurements from e"*" e~ collisions and lattice 
simulation results, respectively. Our analysis of the theoretical uncertainties is based on 
different implementations of the perturbative series and on independent variations of the 
renormalization scales for the mass and the strong coupling. Taking into account the 
resulting set of series to estimate perturbative uncertainties is crucial, since some ways 
to treat the perturbative expansion can exhibit extraordinarily small scale dependence 
when the two scales are set equal. As an additional refinement, we address the issue that 
double scale variation could overestimate the perturbative uncertainties. We supplement 
the analysis with a test that quantifies the convergence rate of each perturbative series by 
a single number. We find that this convergence test allows to determine an overall and 
average convergence rate that is characteristic for the series expansions of each moment, 
and to discard those series for which the convergence rate is significantly worse. We obtain 
rficijnc) = 1.288 ± 0.020 GeV from the vector correlator. The method is also applied to 
the extraction of the MS bottom quark mass from the vector correlator. We compute the 
experimental moments including a modeling uncertainty associated to the continuum region 
where no data is available. We obtain mb{rnb) = 4.176 ± 0.023 GeV. 



Contents 


1 Introduction 2 

2 Theoretical Input 5 

2.1 Perturbative Contribution 5 

2.2 Gluon Condensate Contribution 8 

2.3 Ratios of Moments 9 

3 Previous Results and Scale Variations 11 

4 Convergence Test 14 

5 Lattice Simulation Data 18 

6 Computation of the Experimental Moments for the Bottom Correlator 19 

6.1 Narrow Resonances 20 

6.2 Threshold Region 21 

6.2.1 Subtraction of the T(4S') Radiative Tail 21 

6.2.2 Deconvolution of Initial-State Radiation 22 

6.2.3 Determination of the Unfolding Error Matrix 23 

6.2.4 Contribution of the Threshold Region 24 

6.3 Continuum Region 25 

6.4 Final Results for the Experimental Moments 27 

7 Comparison to other Determinations of the Experimental Moments 28 

8 Computation of the Experimental Values for the Ratios of Moments 29 

9 Results 30 

9.1 Results for the Charm Mass from the Vector Correlator 32 

9.2 Results for the Charm Mass from the Pseudoscalar Correlator 34 

9.3 Results for the Bottom Mass from the Vector Correlator 36 

10 Comparison to other Determinations 38 

10.1 Charm Mass 39 

10.2 Bottom Mass 40 

11 Conclusions 41 

A Numerical Values for the Perturbative Coefficients 43 


- 1 - 




1 Introduction 


Precise and reliable determinations of the charm and bottom quark masses are an important 
input for a number of theoretical predictions, such as Higgs branching ratios to charm and 
bottom quarks or for the corresponding Yukawa couplings [1, 2], They also affect the 
theoretical predictions of radiative and inclusive B decays, as well as rare kaon decays. For 
example, the inclusive semileptonic decay rate of B mesons depends on the fifth power of 
the bottom quark mass. These weak decays provide crucial methods to determine elements 
of the CKM matrix, which in turn are important for testing the validity of the Standard 
Model, as well as for indirect searches of new physics. In this context, having a reliable 
estimate of uncertainties for the quark masses is as important as knowing their precise 
values [3]. Due to confinement quark masses are not physical observables. Rather, they are 
scheme-dependent parameters of the QCD Lagrangian which have to be determined from 
quantities that strongly depend on them. 

One of the most precise tools to determine the charm and bottom quark masses 
is the QCD sum rule method, where weighted averages of the normalized cross section 
Re+e-^gq+X, with q = C,b, 


~ / „n+l e-^ qq +X [S) , Rg+e-^ qq +X [S) — , ("si ’ 

J ^ ^e+e“^ 

can be related to moments of the quark vector current correlator liy [4, 5]: 


( 1 . 1 ) 


= 


rV,th _ d" 


n! 


ds*^ 


ny(s) 


s=0 


f{x) = qix)'yf^q{x ), 


{9^luS- q^qu)llv{s) = ~^ J da; e*'?"' ( 0 |r j^(x)ji.(0)| 0; 


( 1 . 2 ) 


Here Qq is the quark electric charge and ^/s = ^de e'’'e“ center-of-mass en¬ 

ergy. Given that the integration over the experimental ii-ratio extends from the quark 
pair threshold up to infinity but experimental measurements only exist for energies up to 
around llGeV, one relies on using theory input for energies above that scale (which we 
call the “continuum” region). For the charm moments, the combination of all available 
measurements is actually sufficient to render the experimental moments essentially inde¬ 
pendent of uncertainties one may assign to the theory input for the continuum region [6]. 
For the bottom moments, the dependence on the continuum theory input is very large, 
and the dependence of the low-n experimental moments on unavoidable assumptions about 
the continuum uncertainty can be the most important component of the error budget, see 
e.g. [7]. In fact, the use of the first moment MY to determine the bottom mass appears to 
be excluded until more experimental data becomes available for higher energies. 

Alternatively one can also consider moments of the pseudoscalar current correlator to 
extract the heavy quark masses. Experimental information on the pseudoscalar correlator 
Hp is not available in a form useful for quark mass determinations, but for the charm quark 
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very precise lattice calculations have become available recently [8]. For the pseudoscalar 
correlator it turns out that the first two Taylor coefficients in the small-g^ expansion need 
to be regularized and defined in a given scheme, and that the third term (which we will 
denote by Mq) is hardly sensitive to niq. We adopt the definitions 


Ilp{s) = i / dxe*''^(0|rjp(x)jp(0)|0) , jp{x) = 2mqi q{x)'y 5 q(x), 


(1.3) 


127r^Q^ d 


^p,th ^ i—P{s) 

n! ds” ^ ^ 


s=0 


P{s) = 


np(s) ~ np(o) — np(o) s 


where the explicit mass factor in the definition of the pseudoscalar current ensures it 
is renormalization-scheme independent. 

For small values of n such that mq/n > Aqcd; the theoretical moments for the vector 
and pseudoscalar correlators can be computed in the framework of the OPE, i.e. as an 
expansion in vacuum matrix elements involving operators of increasing dimension [4, 5]. 
The leading matrix element corresponds to the perturbative QCD computations, which 
greatly dominates the series. Nonperturbative corrections are parametrized by vacuum 
condensates, and we find that even the leading correction, given by the gluon condensate, 
has a very small effect for low n, particularly so for the bottom correlator. For moments 
at low values of n, it is mandatory to employ a short-distance mass scheme such as MS [9], 
which renders the quark mass fnq{^rn) dependent on its renormalization scale similar 
to the strong coupling as{^a)i which depends on Ha- This method of determining the heavy 
quark masses with high precision is frequently called relativistic charmonium/bottomonium 
sum rules. 

For the perturbative term, the exact analytic expressions for the II functions are known 
at 0{ag) and 0{as), [10]. Therefore any moment can be obtained simply by Taylor ex¬ 
panding around q"^ = 0. At 0{ag) moments are known to up to ra = 30 [11-15]. At 
0{ag) they are known analytically for n = 1 [16-18], n = 2 and n = 3 (and even n = 4 
for the pseudoscalar correlator) [19, 20]. Higher moments at 0{ag) have been determined 
by a semi-analytical procedure [21-24]. The Wilson coefficient of the gluon condensate 
contribution is known to 0(as) [25]. 

The most recent determinations of the MS charm mass from charmonium sum rules 
for the vector correlator [6, 26, 27] obtain very accurate results, but differ in the way 
they estimate theoretical uncertainties, and also in the computation of the moments from 
experimental data. Concerning the estimate of the perturbative uncertainties. Ref. [6] 
obtained 19MeV compared to 1 and 2MeV obtained in Refs. [26] and [27], respectively. 
The discrepancy arises from two differences. First, in Refs. [26, 27] the renormalization 
scales /Xrm and /Tq were set equal, while in Ref. [6] it was argued that they should be 
varied independently. Second, in Refs. [26, 27] the lowest renormalization scale was chosen 
to be 2GeV, while in Ref. [6] variations down to the charm mass value were considered. 
For the case of the pseudoscalar moments, the HPQCD collaboration has made a number 
of very accurate predictions for charm and bottom masses [8, 28-30], the last of which 
has the smallest uncertainty claimed so far for the charm mass. In all these analyses the 
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renormalization scales and //q, are set equal when estimating the truncation uncertainty. 
A detailed discussion on the estimates of theoretical and experimental uncertainties can be 
found in Secs. 3 and 10 of this article (see also Ref. [6]). 

Similarly, bottomonium sum rules have been used to determine the bottom mass from 
low-n moments. To the best of our knowledge, the most recent and precise determina¬ 
tions are from Refs. [27, 31]. These two analyses estimate their perturbative uncertainties 
in the same way as the corresponding charm mass extractions from the same collabora¬ 
tions [26, 27]. Furthermore, when it comes to compute the experimental moments, they use 
theoretical input at 0{a^) with perturbative uncertainties to model the high-energy region 
(continuum region) of the spectrum. As we discuss in this work, similar caveats as for their 
charm analyses can be argued to also affect their bottom quark results. 

In this work we revisit the charmonium sum rules for the vector correlator, refining 
our perturbative error estimate from Ref. [6] by incorporating a convergence test. The 
convergence test addresses the issue that the independent variation of and /Xq, together 
with the relatively large value of the Og close to the charm mass scale, might lead to an 
overestimate of the perturbative uncertainty. The convergence test allows to quantify the 
convergence property of each perturbative series with a single parameter and to discard 
series for which the convergence is substantially worse then for the rest of the series. We 
show that this procedure is meaningful, since the complete set of series for the moments 
shows a strongly peaked distribution in these convergence values, which allows to define 
an overall convergence for the set of perturbative series. This leads to a reduction of 
the perturbative uncertainty from 19 to 14MeV, and the corresponding result for the MS 
charm mass supersedes the main result given in Ref. [6]. We also apply this improved 
method of estimating theory uncertainties to obtain a new charm mass determination from 
the pseudoscalar correlator, and to extract the bottom mass from the vector correlator. For 
the latter, we compute the bottom experimental moments by combining contributions from 
narrow resonances, experimental data taken in the continuum, and a theoretical model for 
the continuum region. We carefully study the assignment of adequate uncertainties to this 
last contribution, to make sure that the the model dependence is reduced to an acceptable 
level. 

This paper is organized as follows: In Sec. 2 we summarize the theoretical framework 
introduced in [6], and adapt it to cover the case of the pseudoscalar moments. We also 
introduce the ratios of moments, also used before in Ref. [9], and the perturbative expansions 
associated to them. Sec. 3 contains a brief summary of the results obtained in [6], and the 
discussion is extended to the case of the pseudoscalar correlator and the bottom mass. In 
Sec. 4 we introduce the convergence test, and discuss how it allows to identify and discard 
series with a bad convergence. Sec. 5 contains a discussion on the lattice simulation results 
we use for our analysis. In Sec. 6 we present our computation of the experimental moments 
for the bottom correlator. The results are compared to previous determinations in Sec. 7. 
The computation of the ratio of experimental moments is presented in Sec. 8. Our final 
results for the quark masses are given in Sec. 9. The results are compared to previous charm 
and bottom mass analyses in Sec. 10. We present our conclusions in Sec. 11. In Appendix A 
the numerical values of the perturbative coefficients that enter into our analysis and are 
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not yet provided by Ref. [6] are collected for the convenience of the reader. 

2 Theoretical Input 

2.1 Perturbative Contribution 

The moments of the vector and pseudoscalar current correlators are defined in Eqs. (1.1) 
and (1.3), respectively. In the OPE framework they are dominated by the perturbative 
contribution (that is, a partonic computation), which exhibits a nonlinear dependence on 
the mass. Within perturbation theory one can decide to manipulate the series expansion 
to get a more linear dependence on the mass. Conceptually there is no preference. As 
advocated in our previous analyses [6] , one might consider various versions of the expansion 
to reliably estimate the perturbative uncertainties. Four types of expansion were suggested 
in Ref. [6], which we briefly review below. 

(a) Standard fixed-order expansion 

We write the perturbative vacuum polarization function as 


^x{s,nf,CVs (l^o)) (^m)) A^m) — 


1 


oo 




X 


( 2 . 1 ) 


where X = V, P for vector and pseudoscalar currents, respectively. Note that for 
notation reasons, in Eqs. (2.1, 2.5, 2.7, 2.8) we use Ylp{(f‘) = P{q‘^), where P is the twice- 
subtracted pseudoscalar correlator defined in Eq. (1.3). Here Qq is the quark electric charge 
with q = c^b, and n/ = 4,5 for charm and bottom, respectively. In full generality, the 
perturbative moments can be expressed as the following sum: 


= 


[4mPHm)y 
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z,a,6 


as (a^o 
vr 


[Cx{nf)t’’’ ln“ 


Jn,2 


mq 


il^m) 




ln'= 


mq 






.( 2 . 2 ) 


This is the standard fixed-order expression for the perturbative moments. The numer¬ 
ical values for the [Cv{nf = 5)]“’^ coefficients are given in Table 8 (for the vector current 
with Uf = 4, these coefficients can be found in Table 1 of Ref. [6]). Likewise, [Cp{nf = 4)]“’^ 
are collected in a numerical form in Table 11. 

The expression in Eq. (2.2) is the common way to write the perturbative series of the 
moments. However, as noted in Ref. [6], the nonlinear dependence on friq of the standard 
fixed-order expansion has the disadvantage that for charm (bottom) quarks there are fre¬ 
quently no solutions for the mass in the sum rule mass determination, for moments higher 
than the first (second), for some set of values of the renormalization scales. 


(b) Linearized expansion 
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One can linearize the the fixed-order form expansion of Eq. (2.2) with respect to the 
exponent of the quark mass pre-factor by taking the 2n-th root. This choice is e.g. made 
in Ref. [28], and in general one can write: 



l/2n. 


1 


Ih 




(2.3) 


The coefficients \Cv{nf = 5)]“’^ and [Cp{nf = 4)]“’^ are given in Tables 9 and 12, 
respectively (for n/ = 4 the coefficients for the vector current can be found in Table 2 
of Ref. [6]). Even though relation (2.3) still exhibits some nonlinear dependence on friq 
through perturbative logarithms, we find that it always has a solution for the quark mass. 

(c) Iterative linearized expansion 

For the expansion methods (a) and (b) shown in Eqs. (2.2) and (2.3), one solves for 
the quark masses mc,b(/^m) numerically keeping the exact mass dependence on the theory 
side of the equation. Alternatively, one can solve for iteratively order by order, 

which is perturbatively equivalent to the exact numerical solution, but gives different nu¬ 
merical results. The method consists of inserting the lower order values for mcfiipm) in 
the higher order perturbative coefficients, and re-expanding consistently. This method has 
been explained in detail in Sec. 2.1(c) of Ref. [6], and we only quote the final results here: 

where the numerical value of the coefficients [Cyi'iT'f = ^)]n\ [Cp{nf = 4)]“’^ are 
collected in Tables 10 and 13, and the values for the vector current with n/ = 4 can be 
found in Table 3 of Ref. [6]. 

By construction, the iterative expansion always has a solution for the quark mass. 
Accordingly, potential biases on the numerical analysis related to any possible nonlinear 
dependence are eliminated. 

(d) Contour-improved expansion 

For the expansions (a), (b) and (c) the moments and the quark masses are computed for 
a fixed value of the renormalization scale Using the analytic properties of the vacuum 
polarization function, one can rewrite the fixed-order moments as integrals in the complex 
plane. This opens the possibility of making dependent on the integration variable. 


i,a^b 

1 




TT 


In' 


—{ 0)2 

rriq 

Pm 


In'^ 


pI 


(2.4) 
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Figure 1. One possible integration path in the complex s-plane for the computation of the 
contour-improved moments. 


in analogy to the contour-improved methods used for r-decays (see e.g. Refs. [32-37]). 
Therefore we define the contour-improved moments [38] as (see Fig. 1), 


mx,c ^ 


ds 


nx(s, Uf, nig)), mq{nm), rng),iJ.m), (2.5) 


™2\\ — 




Ic ■S' 


,n+l 


and we employ the following path-dependent first used in Ref. [38] 


= ( 2 - 6 ) 

It weights in a different way the threshold versus the high energy parts of the spectrum. 
It was shown in Ref. [6] that the resulting moments Mn ’ can be obtained analytically 
from the Taylor expansion around s = 0 of the vacuum polarization function using an 
s-dependent /i^(s,mg): 


OO 

^X^[s,a^r^\filis,m^)),mq{flrn),fia{s,fn^),flm) = ^ s" . (2.7) 

This trick works because as{iJ.a{s,nlg )) has the same cut as the fixed-order expression 
for fix- Other choices could spoil this property. Expanding the analytic expression for 
Mn ’ on as at a given finite order, one recovers the fixed-order moments . This shows 
that the dependence on the contour is only residual and represents an effect from higher 
order terms from beyond the order one employs for the calculation. 

The contour-improved moments have a residual sensitivity to the value of the vacuum 
polarization function at zero momentum transfer.^ For the case of the vector correlator this 
value depends on the UV-subtraction scheme and corresponds to 11(0) = For the case 
of the pseudoscalar correlator, Mq is scheme-independent, since P{q^) already includes two 
UV subtractions. However, one could as well define a three-times-subtracted pseudoscalar 

^This means that the dependence vanishes in the large-order limit. 
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correlator of the form P{(f‘) = Piq^) —-P(O). Slightly abusing notation, we denote P as the 
“on-shell” scheme for P{q‘^), and the twice subtracted (original) definition as the MS scheme 
for P{q‘^). Using the OS scheme with 11^(0) = 0 for either vector or pseudoscalar correlator, 
we find that the first moment for the contour-improved expansion gives exactly the first 
fixed-order moment, . Thus, in order to implement a non-trivial modification, 

and following Ref. [6], we employ the MS scheme for 11^(0) defined for /i = fnq{rnq), and 
the twice-subtracted expression for P{(f‘). Generically it can be written as 


nr(0,n/) = 


E 

i,a,b 


TT 


[Cx{nf)roi ln“ 


_ 9 


(Pn 




In* 


_ 9 

mf: 


{pn 


Pc 


( 2 . 8 ) 


The numerical values for the coefficients are collected in Table 7 for the vector 

correlator with 5 flavors and the pseudoscalar correlator with 4 flavors. In Table 4 or Ref. [6] 
one finds the the numerical values of [C'y(nj = d)]^’!*. 


2.2 Gluon Condensate Contribution 

We estimate nonperturbative power corrections by including the gluon condensate contri¬ 
bution. The gluon condensate is a dimension-4 matrix element and gives the leading power 
correction in the OPE for the moments [39, 40] 


+ ... (2.9) 

Here the ellipses represent higher-order power corrections of the OPE involving con¬ 
densates with dimensions bigger than 4. The Wilson coefficients of the gluon condensate 
corrections are known to 0(as} accuracy [25]. Following Ref. [41], we express the Wilson 
coefficient of the gluon condensate in terms of the pole mass, since in this way the correction 
is numerically more stable for higher moments. However, as we did in Ref. [6], we still write 
the pole mass in terms of the MS quark mass at one loop. The resulting expression reads 




Mq 


1 

(4M2)^+2 \ tt / RGI 




rriqifim) 



a^r^\pa) 

TT 



( 2 . 10 ) 


We use the renormalization group invariant (RGI) scheme for the gluon condensate 
[42]. The numerical value of the [ay(n/ = 5)]“ and [ap{nf = 4)]“ coefficients are collected 
in Table 6. The values for [ay(nj = 4)]“ can be found in Table 5 of Ref. [6]. For methods 
(b) and (c) one can obtain the gluon condensate contribution by performing simple algebra 
operations and re-expansions in and (G^). For method (d) we employ Eqs. (2.9) and 
(2.10) as shown. For the RGI gluon condensate we adopt [43] 


^^2 


= 0.006 ± 0.012 GeV"^. 
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PiYf^rl r\rr\(^r 


Linearized 




(a) 


(b) 




(c) (d) 


Figure 2. Contour plots for fnc{mc) as obtained from the first moment of the pseudoscalar cor¬ 
relator Ml , as a function of /i^ and Hm at 0{a^), for methods (a)-(d). The shaded areas rep¬ 
resent regions with < tric(TOc), and are excluded of our analysis. For this plot we employ 

oisijnz) = 0.118. 


2.3 Ratios of Moments 

An alternative set of observables which are also highly sensitive to the quark masses are 
the ratios of consecutive moments. To that end we define R^{nf) = M^_^i {nf)/M^{nf). 
Such ratios are proportional to the inverse square of the quark mass for any value of n. 
Their perturbative series can be expressed as an expansion in powers of analogous to 
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Fixed order Linearized 



(a) 

Iterative 


(b) 

Contour Improved 



(c) 


(d) 


Figure 3. Contour plots for fnb{fnb) as obtained from the second moment of the vector correlator 
M 2 with n/ = 5, as a function of and fim at ©(af), for methods (a)-(d). The shaded areas 
represent regions with /Tq < fnb{mb), and are excluded of our analysis. For this plot we employ 
as{mz) = 0.118. 


Eq. (2.2), with the replacements [4m^(//m)]" -t (/tm) and [^^(n/)]^’^^ 

Their computation is trivial, as one only needs to take the ratio of the two consecutive 
theoretical moments and re-expand as a series in We call this the standard fixed- 

order expansion, analogous to the expansion (a) of Sec. 2.1. The numerical expressions for 
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the [Rx{nf)f'^\ coefficients for the vector correlator with n/ = 4,5 are given in Table 14, 
and for the pseudoscalar correlator with n/ = 4 in Table 15. As for the regular moments, 
we find that the nonlinear dependence of on the quark mass sometimes causes that 
there is no numerical solution for nlq. 

By taking the square root of the ratio of two consecutive moments one gets a linear 
dependence on the inverse of the quark mass. The corresponding theoretical expression is 
obtained by re-expanding the perturbative expansion of (uf) as a series in powers 

of Thus we obtain an expression of the form of Eq. (2.3) with the replacement 

[Cx(n —7- [Rx{nf)y^^fj. This is referred to as the linearized expansion, in analogy to 
the expansion (b) of Sec. 2.1. The numerical values for the [Rxinf)Y^\ coefficients are 
collected for the vector correlator with n/ = 4,5 in Table 16, and for the pseudoscalar 
correlator with n/ = 4 in Table 17. 

Finally, one can use y^R^{nf) to solve for mq{iirn) in an iterative way, exactly as 
explained in Sec. 2.1 for expansion (c). The theoretical expression is analogous to Eq. (2.4) 
with the replacement [Cx{nf)f'^\ —>• [Rxinf)]]^\- We collect the numerical values for the 
[Rxinf)]]fij coefficients, in Tabs. 18 and 19 for the vector {rif = 4,5) and pseudoscalar 
{uf = 4) correlators, respectively. We call this the iterative linearized expansion. 

One cannot implement a contour-improved expression for the ratios of moments, as 
they cannot be computed as the contour integral of a correlator. For the ratios of moments, 
in any of the three expansions, one can include non-perturbative corrections in the form of a 
gluon condensate OPE term, just using Eq. (2.10) and performing simple algebra operations 
and re-expansions in and (G^). 

3 Previous Results and Scale Variations 

In a number of recent low-n sum rule analyses [8, 17, 26-29, 31, 44], which determined 
the charm and bottom quark masses with very small uncertainties using 0(af) theoretical 
computations for the moments [21-24], the theory uncertainties from the truncation of 
the perturbative series have been estimated with the scale setting firn = fJ-a based on 
just one type of expansion, which was either the fixed-order [expansion (a)] for the vector 
correlator, or the linearized [expansion (b)] for the pseudoscalar correlator. In Ref. [6] we 
analyzed the perturbative series for the moments MY 2 3 4 of th® charm vector correlator 
at 0{al) using an alternative way to estimate the perturbative uncertainties, based on 
the four different expansion methods (a)-(d), as explained in Sec. 2.1. We also focused 
on the question whether renormalization scale variation restricted to = fJ-a leads to a 
compatible estimate of the perturbative uncertainties. From our analysis we found: 

• The extractions for the MS charm mass using the expansions (a) - (d) with correlated 
variations of Hm and /Tq, (e.g. = fJ^a) for the vector correlator can lead to very small 

scale variations, which can be very different depending on the method.^ Moreover, 

^We judge the compatibility of the perturbative error estimates based on the size of scale variations 
alone, i.e. without accounting at the same time for other sources of uncertainties such as experimental 
errors or the uncertainty in the value of the strong coupling. 
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for some expansions also the results from the different orders can be incompatible to 
each other. It was therefore concluded that using correlated scale variation and one 
type of expansion can lead to an underestimate of the perturbative uncertainty. 

• Uncorrelated (i.e. independent) variation of and leads to charm mass extrac¬ 
tions with perturbative uncertainty estimates that are in general larger, but fully 
compatible among the expansions (a) - (d) and for the different orders. It was there¬ 
fore concluded that fj,m and should be varied independently to obtain a reliable 
estimate of the perturbative uncertainty. 

• The size of the charm mass perturbative uncertainty has a significant dependence on 
the value of the lower bound of the range of the scale variation. The choice of the 
upper bound has a marginal impact. 

• The pattern of size of the correlated scale variations for the different expansions can 
be traced back to the form of the contours of constant charm mass in the 
plane, which happens to be located along the diagonal ~ ^^a for expansions (a) 
and (b), but roughly orthogonal for expansions (c) and (d), see Fig. 6 of Ref. [6]. 

For example,^ in Ref. [27] method (a) has been used for MY with = /^q varied 
between 2 and 4 GeV, quoting a perturbative error estimate of 2 MeV. For the expansion 
methods^ (a) - (d) we obtain for the same scale variation 1.2893 ± 0.0007, 1.2904 ± 0.0004, 
1.2963 ± 0.0045 and 1.3009 ± 0.0020 GeV, respectively, for rndmc), which are inconsistent. 
This can be compared to the corresponding results using independent variations as sug¬ 
gested in Ref. [6]. Using 2 GeV < IJ^m < 4 GeV we obtain 1.291 ± 0.003, 1.291 ± 0.003, 
1.296 ± 0.005 and 1.302 ± 0.003GeV, respectively, for expansions (a)-(d). These results 
are not consistent either. It was furthermore argued in Ref. [6] that an adequate variation 
range should include the charm mass itself (after all, that is the scale that governs the 
series), motivated by the range 2me ± rric around the pair production threshold. Thus, 
adopting independent scale variation in the range rueijne) < fJ-m, < 4 GeV one obtains 
1.287 ± 0.018, 1.287 ± 0.015, 1.282 ± 0.019 and 1.291 ± 0.014 GeV respectively. The re¬ 
sults show consistency and demonstrate the strong dependence on the lower bound of the 
renormalization scale variation. The outcome is illustrated graphically in Fig. 4(a), and 
the order-by-order dependence in Fig. 1 of Ref. [45]. In Ref. [6] we also explored scale 
setting in which fim was fixed to mcijfic) and only /Tq was varied. The outcome is shown 
in Figs. 4 and 5 of that reference. The contour lines in the plane for the mass 

extraction from the first moment of the vector correlator for all methods are shown in 
Fig. 6 of Ref. [6]. The final result quoted in Ref. [6], using as{mz) = 0.1184 ± 0.0021, was 

®The size of the scale variations quoted in this paragraph applies to fficime) as well as to mc(3GeV), 
and all numerical results are obtained at ©(a®). We also stress that there are no perturbative instabilities 
concerning the use of the RGE down to the scale mc(mc). 

^To compute the charm and bottom masses in this section we use = 0.2121 GeV“^ for the charm 

vector correlator, result obtained in Ref. [6], = 0.1402GeV“^ for the charm pseudoscalar correlator, 

from Ref. [8], and our own computation = 2.834 x 10“® GeV”"^ for the bottom vector correlator, 

see Sec. 6. We also use as{mz) = 0.1184. 
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Figure 4. Charm and bottom mass values from the first [second] moment of the vector (a) for 
charm [(c) for bottom] and pseudoscalar [(b), charm] currents at ©(af); and for the ratio of the 
second over the first moment for the vector [(d) for charm, (g) for bottom] and pseudoscalar [(e), 
charm] correlators. We show the outcome of various scale variations for the perturbative expansions 
(a)-(d) [(a) “(c) for ratios], where green (rightmost) corresponds to 2 GeV < = Ma < 4GeV 

]5GeV < fim = < 15 GeV for bottom], blue (second from the left) 2 GeV < Hm, fJ-a < 4 GeV 
[5 GeV < ^J-a < 15 GeV for bottom], purple (second from the right) mcirric) < fJ-a < 4GeV 
^mbirfib) < MmiMa < 15 GeV for bottom] and in red (rightmost) we supplement the latter variation 
with a cut on the series with larger values of V),- 


mc{mc) = 1.282 ± (0.006)stat ± (0.009),yst ± (0.019)pert ± (0.010)„, ± (0.002 )(gg> GeV 
based on the iterative expansion method (c). 

We have repeated this analysis for the first moment of the pseudoscalar correlator 
Mf. Ref. [8] uses method (b) with the same scale variation as Ref. [27], quoting 4MeV 
for the truncation error. For methods (a) - (d) and using 2 GeV < fj,m = l^a ^ ^ GeV 
we obtain 1.276 ± 0.003, 1.277 ± 0.004, 1.275 ± 0.005 and 1.297 ± 0.004 GeV, respectively. 
For independent double scale variation between 2 and 4 GeV we obtain, 1.276 ± 0.013, 
1.277±0.012, 1.271±0.012 and 1.294±0.012 GeV, and if we use mc(mc) as the lower bound 
to we obtain 1.260±0.039, 1.267±0.037GeV, 1.259±0.041 and 1.272±0.034. These results 
are displayed graphically in Fig. 4(b). The contour lines for the mass extraction from the 
first moment of the pseudoscalar correlator for all methods are shown in Fig. 2. We see 
that the results show a qualitative agreement with the situation for the vector current, but 
at a level of perturbative scale variations that are in general roughly larger by a factor of 
two. 

A similar study can be performed for the extraction of the bottom mass rnb{fnb) from 
the second moment of the vector correlator ■ Ref. [27] uses the fixed-order expansion 
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[method (a)] and correlated scale variation between 5 GeV < < 15 GeV, quoting 

a perturbative error of 3MeV. For the same variation we obtain 4.1781 ± 0.0005, 4.1771 ± 
0.0015, 4.1818±0.0034 and 4.1792±0.0044GeV for methods (a) and (d), respectively. As in 
the charm case the results are not consistent, but the variations of the results have a much 
smaller size, as is expected from the fact that for the bottom the renormalization scales 
are much larger. For independent variation between the same values we get 4.183 ± 0.008, 
4.181 ± 0.006, 4.180 ± 0.006 and 4.186 ± 0.013 GeV. Finally, if the lower limit of the double 
variation starts at we find 4.179 ± 0.011, 4.181 ± 0.011, 4.175 ± 0.011 and 4.184 ± 

0.0015 GeV. These results are collected in Fig. 4(c). The corresponding fnb{fnb) contours 
in the plane are shown in Fig. 3. As for the charm case, we find fully consistent 

results for the independent scale variation and using as the lower bound. 

We have also studied the ratio of the first over the second moments for the three cases, 
and observe a very similar pattern. We do not provide a detailed discussion in the text, 
but display the outcome graphically in Figs. 4(d) to 4(f). 

4 Convergence Test 

At this point it is useful to consider the perturbative series for all choices of /Xq, and 
as different perturbative expansions, which can have different convergence properties. To 
estimate the perturbative uncertainties one analyzes the outcome of this set of (truncated) 
series. While the uncorrelated scale variation certainly is a conservative method, one pos¬ 
sible concern is that it might lead to an overestimate of the size of the perturbative error. 
For instance, this might arise for a non-vanishing value of ln(/Xm/A*a) in connection with 
sizeable values of as{na) foi' fJ-a close to the charm mass scale, which might artificially 
spoil the convergence of the expansion. One possible resolution might be to simply reduce 
the range of scale variation (such as increasing the lower bound). However, this does not 
resolve the issue, since the resulting smaller variation merely represents a matter of choice. 
Furthermore, there is in general no guarantee that the series which are left have a better con¬ 
vergence despite the fact that the overall scale variation might become reduced. Preferably, 
the issue should be fixed from inherent properties of the perturbative series themselves. It 
is possible to address this issue by supplementing the uncorrelated scale variation method 
with a convergence test constraint, which we explain in the following. 

We implement a finite-order version of the root convergence test. Let us recall that in 
mathematics, the root test (also known as Cauchy’s radical test) states that for a series of 
terms a„, S[a] = Yln^n, if the quantity Voo defined as 

Hoc = limsup(an)^/" , (4.1) 

n^oo 

is smaller (bigger) than 1, the sum is absolutely convergent (divergent). If 14o ap¬ 
proaches 1 from above then the series is still divergent, otherwise the test is not conclusive. 
In Eq. (4.1) limsup stands for the superior limit, which essentially means that in case of 
oscillating series, one takes the maximum value of the oscillation. In the context of our 
analysis with truncated series the relevant property is that a smaller l^o implies a better 
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convergent series. For the different expansion methods we use, it is simplest to apply the 
method directly to the sequence of quark masses that are extracted order by order, rewrit¬ 
ing the results as a series expansion. Since we only know a finite number of coefficients of 
the perturbative series, we need to adapt the test. We now introduce I 4 and proceed as 
follows: ^ 

(a) For each pair of renormalization scales (/imj /^o) we define the convergence parameter 
Vc from the charm mass series rricijnc) = -|- 5rrS^'^ + resulting from 

the extractions at 


rjmd) 

c = j ’ (. ^(0) ) 

(b) The resulting distribution for Vc values can be conveniently cast as a histogram, and 
the resulting distribution is a measure for the overall convergence of the perturba¬ 
tive expansion being employed. We apply the convergence analysis to the region 
VLcirfic) < < 4GeV for charm, and mf,{fnh) < Ha, IJ^m < 15GeV for bottom. If 

the distribution is peaked around the average (I 4 ) it has a well-defined convergence. 
Hence discarding series with Vc ^ (Vc) (particularly if they significantly enlarge the 
estimate of the perturbative error) is justified. 


Fig. 5(a) shows the Vc distributions for expansions (a)-(d) for the extraction of the 
charm mass from the vector moment MY ■ We find (fc)double = (0.15, 0.15, 0.17, 0.19)® and 
that the distributions are peaked around {Vc), indicating a very good overall convergence. 
The scale variation error (defined as half the overall variation) as a function of the fraction 
of the series (with the largest I 4 values) that are being discarded is shown in Fig. 6(a). We 
see that only around 2% of the series with the highest Vc values by themselves cause the 
increase of the scale variation from well below 15MeV to up to 20MeV. These series are 
located at the upper-left and lower-right corners of Figs. 2 and 3, and Fig. 6 of Ref. [6], 
corresponding to values of Hm and fj-a far from each other. Given that these series have 
very large Vc values and do not reflect the overall good convergence behavior of the bulk of 
the series, it is justified to remove them from the analysis. 

The Vc distributions for the pseudoscalar first moment Mf are shown in Fig. 5(b), again 
showing a clear peak. However, with (V);)doubie = (0.24, 0.24, 0.25, 0.21) [for correlated 
variation (14)corr = (0.22, 0.23, 0.22, 0.15) with 2GeV as the lower bound], the average 
Vc values are significantly larger than for the vector correlator, indicating that the overall 
perturbative convergence for the pseudoscalar moment is still excellent but worse than for 
the vector moment. This means that the vector correlator method is superior, and we 

®One could think of implementing the ratio test as well. However, since we only known a small number 
of terms, it is likely that one of them becomes close to zero, making one of the ratios blow up. This makes 
this test very unstable. 

^Interestingly, the same analysis for the correlated variation with 2 GeV < fJ-a = Hm < 4 GeV yields 
(Tcjcorr = (0.14, 0.16, 0.19, 0.19), which is similar to the outcome for the double variation. This same 
observation can be made for the rest of the correlators. 
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Figure 5. distribution for fficimc) from the first moment of the vector (a) and pseudoscalar 
(b) correlator, and for fFib{mb) for the second moment of the vector correlator (c), for expansions 
(a)-(d). The three lower panels show the same for the ratio of the second over the first moment 
for expansions (a)-(c). 


expect that the perturbative uncertainty in the charm mass from the pseudoscalar is larger. 
This expectation is indeed confirmed as we discussed in Sec. 3, see also Sec. 9. Fig. 6(b) 
shows that the effect of discarding the series with the worst convergence is very similar to 
that of the vector correlator. 

For our determination of the bottom mass we use the second moment (see Sec. 6 for 
a discussion on why we discard the first moment), and employ uncorrelated scale variations 
in the range mbinib) < IJ-m-, IJ^a < 15 GeV. Fig. 5(c) shows the corresponding histograms, and 
we find that the convergence test yields (V);)double = (0.13, 0.11, 0.12, 0.15) for expansions 
(a) - (d) [for the correlated variation with scales set equal and 5 GeV < Ha = IJ-m < 15 GeV 
we find (V[;)corr = (0.13, 0.09, 0.13, 0.15)]. As expected, the averages for the bottom are 
much smaller than for the charm. We further find that discarding series with the highest 
Vc values only has minor effects on the perturbative error estimate for fractions up to 5%, 
see Fig 6(c). This is a confirmation that the series for bottom moments overall are more 
stable, which is again expected from the fact that perturbation theory should work better 
for the bottom than for the lighter charm. 

The behavior of the ratios of moments is very similar as that for regular moments, as 
can be seen in Figs. 6 and 5, panels (d) - (f). We find the following average values for Vc for 
methods (a) - (d): ratios of charm vector moments (V'c)double = (0.19, 0.18, 0.19) [(Vc)corr = 
(0.16, 0.16, 0.23)]; ratios of charm pseudoscalar moments (V'c)doubie = (0.25, 0.23, 0.18) 
[(lc)corr = (0.25, 0.20, 0.16)]; ratios of bottom vector moments (V'c)doubie = (0.13, 0.12, 0.13) 
[(lc)corr = (0.13, 0.11, 0.14)]. Therefore we conclude that the perturbative convergence of 
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Figure 6. Half of the scale variation of mq{rnq) at 0{al) as a function of the fraction of the 
discarded series with highest 14 values for the first moment of the vector (a) and pseudoscalar (b) 
correlators for charm, the second moment of the vector correlator for bottom (c); the ratio of the 
second over the first moment for the vector [charm (d) and bottom (f)] and pseudoscalar [charm 
(e)[ correlators. 


the ratios of moments is in general terms a bit worse than that of regular moments, except 
for the linearized iterative method of the pseudoscalar ratios. 

In our hnal numerical analyses we discard 3% of the series with the worst 14 values. As 
can be seen from Fig. 5, this only affects series with 14 values much larger than the average 
values for the whole set of series. It is our intention to keep the fraction of discarded series 
as small as possible, since it is our aim to remove only series with convergence properties 
that are obviously much worse than those of the bulk of the series. We call this procedure 
trimming in the following. As we see in Figs. 7, the results including the trimming show 
a very good order-by-order convergence for the heavy quark mass determinations, and at 
each order every expansion method gives consistent results for the central values as well as 
for the estimate of the perturbative uncertainties. Figs. 7(a) and 7(b) show the results for 
mc(mc) from the vector and pseudoscalar correlators, respectively, for expansions (a)-(d) 
at and with rndnic) < < 4GeV, using the hrst moment. Figs. 7(d) and 

7(e) show results for methods (a) ~ (c), using the ratio of the second over the hrst moment. 
Analogously, Figs. 7(c) and 7(f) show the results for mft(mb) for the second moment, and 
the ratio of the second over the hrst moment, respectively, with the uncorrelated variation 
mbiruh) < < 15GeV. 
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(a) (b) (c) 



(d) (e) (f) 

Figure 7. Charm and bottom mass values from the first [second] moment of the vector (a) for 
charm [(c) for bottom] and pseudoscalar [(b), charm] currents at and for the ratio of 

the second over the first moment for the vector [(d) for charm, (g) for bottom] and pseudoscalar 
[(e), charm] correlators for expansions (a)-(d) [(a)-(c) for ratios], in red, blue, green and purple, 
respectively. 


5 Lattice Simulation Data 

The pseudoscalar current is not realized in nature in a way which is useful to compute the 
moments of the corresponding correlator from experimental data. Results for the moments 
of the pseudoscalar current correlator can, however, be obtained from simulations on the 
lattice. The strategy of these numerical simulations is to tune the lattice parameters (such 
as bare coupling constant and masses) to a selected number of observables (e.g. the energy 
splitting of T resonances). Once this tuning is performed, the lattice action is fully specified 
and no further changes are implemented. The tuned lattice action can then be used to 
perform all sorts of predictions, moments of correlators among them. Lattice simulations 
have to face a number of challenges, which usually translate into sizeable uncertainties. 
Among those are the extrapolations to the infinite volume and the zero lattice spacing (the 
latter being much harder), as well as the extrapolation to physical quark masses. On top of 
these systematic uncertainties, there are also statistical ones, which are related to the finite 
sampling used to perform the path integrations. On the other hand there are also concerns 
on the type of lattice action which is being used for the fermions. According to Ref. [8], 
the moments of the pseudoscalar correlator are least affected by systematic uncertainties, 
and so HPQCD has focused on those for their subsequent high-precision analyses. 

To the best of our knowledge, HPQCD is the only lattice collaboration which has 
published results on QCD correlators. They have used the staggered-quarks lattice action. 
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and MILC configurations for gluons. These results have been used to determine the charm 
mass and the strong coupling constant [8, 28, 29] with high accuracy, as well as the bottom 
mass [28, 30], with smaller precision. We will use the simulation results as given in Ref. [8], 
even if the results quoted in of [28, 29] are a bit more precise. The reason for this choice is 
that while [8] makes a straightforward extrapolation to the continuum, which is independent 
of the charm mass and ag extractions, in [28, 29] the fit for the quark masses, the strong 
coupling constant and the extrapolation to the continuum is performed all at once, in a 
single fit. Furthermore that fits contains a lot of priors for the parameters one is interested 
in fitting. In any case, as we have seen, the charm mass extraction from the pseudoscalar 
correlator is dominated by perturbative uncertainties, as a result of the bad convergence of 
the series expansion for its moments. 

Ref. [8] provides simulation results for the so-called reduced moments Rk, which are 
collected in their Table II. The index k takes only even values, and starts with the value 
= 4, which is fairly insensitive to the charm mass. Hence the lowest moment we consider is 
Rq. Reduced moments are defined as (up to a global power) the full moment divided by the 
tree-level result. By taking this ratio, the authors of Ref. [8] claim that large cancellations 
between systematic errors take place. The reduced moments are scaleless, and the mass- 
dimension that one obviously needs to determinate the charm quark mass is regained by 
dividing with the mass of the r]c pseudoscalar particle. Thus one can easily translate the 
reduced moments into the more familiar correlator moments with the following relation: 



(5.1) 


where the Cp coefficients correspond to the tree-level terms of the standard fixed-order 
expansion of Eq. (2.2). Although the experimental value for r/c is 2.9836(7) GeV, we use 
the value = 2.980 GeV given in Ref. [8], in order to ease comparison with that analysis. 
In Ref. [29] the value = 2.9863(27) GeV is used. It is claimed that (as for the lattice 
action) it has no QED effects, and the error accounts for cc annihilation. Using the quote 
in Ref. [29] changes Mf by 0.4% and the effect on the charm mass is of the order of 2 MeV. 
The uncertainty in the rjc mass has no effect on the errors. In Table 1 we quote the 
lattice simulation results written as regular moments . 


Mf Aff M[ M[ 


1.402(20) 1.361(40) 1.426(59) 1.558(89) 


Table 1. Lattice simulation results for the moments of the twice-subtracted pseudoscalar correlator 
P{q^) for n/ = 4. Moments given in units of 10“"x GeV“^”. 

6 Computation of the Experimental Moments for the Bottom Correlator 

In this section we present our computation of the moments for the bottom vector current 
correlator from experimental data. These are made of three distinct contributions: the 
narrow resonances below threshold, the region of broader resonances, explored by BABAR 
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[46], and the continuum region, where no data has been taken and some modeling is required. 
The BABAR data has to be corrected for initial-state radiation and vacuum polarization 
effects. In the continuum region we use a model which consists of a combination of a linear 
fit to the BaBar experimental points with energy larger than 11.05 GeV, and perturbation 
theory as a model for missing experimental data, which are joined smoothly by a cubic 
interpolation. We assign a conservative uncertainty guided by the error function of the 
linear fit to the BABAR data in the region with measurements with the highest energy. 
Our determination of the experimental bottom moments differs from Ref. [27] in the way we 
model the uncertainties for the hadronic cross section in the continuum region, plus other 
minor differences in the contributions from the narrow widths and the threshold region, 
see discussion in Sec. 7. We also provide the correlation matrix among different moments, 
which cannot be found in the literature. Therefore, even though there are some similarities 
with the computations outlined in [27], we find it justified to discuss our computation of 
the experimental moments in some detail. 

We note that our results for the moments of the bottom vector current correlator 
have already been used in the analyses of high-n moments of Ref. [47], in the context of 
nonrelativistic large-re sum rules. 


6.1 Narrow Resonances 

The contribution of resonances below the open bottom threshold y/s = 10.62 GeV includes 
T(IS') up to T(45). We use the narrow width approximation to compute their contribution 
to the experimental moments, finding 


= 

iy±n 


OyrTf, 




2n-\-l 


( 6 . 1 ) 


The masses and electronic widths of these four resonances are taken from the PDG 
[48], and the values of the effective electromagnetic coupling constant evaluated at the T 
masses are taken from Ref. [44]. This information is collected in Table 2. We have also 
checked that if one uses a Breit-Wigner instead of the narrow width approximation the 
results change by an amount well within the error due to the uncertainty in the electronic 
width. 

In analogy to what we found in our study of the charm moments [6], the effect of the 
mass uncertainty in the moments is negligible. Therefore one only needs to consider the 
experimental uncertainty in the electronic widths. There is no information on the corre¬ 
lation between the measurements of these widths. The PDG averages of the electronic 
partial widths for the first three resonances is dominated by the GLEO measurement [49].^ 
Therefore we take the approach that half of the width’s uncertainty (in quadrature) is un¬ 
correlated (therefore mainly of statistical origin), whereas the other half is correlated among 
the various resonances (therefore coming from common systematics in the measurements). 


^Refs. [27, 44] assume that the error of the electronic width for the first three narrow resonances is 100% 
correlated, and uncorrelated to that of the T (45). 
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T(15) 

T(25) 

T(35) 

T(45) 

MT(GeV) 

9.46030(26) 

10.02326(31) 

10.3552(5) 

10.5794(12) 

ree(KeV) 

1.340(18) 

0.612(11) 

0.443(8) 

0.272(29) 

/ OQED 1 

Y CK(Af'f ) j 

0.932069 

0.93099 

0.930811 

0.930093 


Table 2. Masses and electronic widths of narrow T resonances [48] and effective electromagnetic 
conpling constant [44]. oqed = 1/137.035999084(51) represents the fine structure constant. 


6.2 Threshold Region 

The region between the open bottom threshold and the experimental measurement of the 
i?fe-ratio at the highest energy, 10.62 GeV < < 11.2062 GeV, is referred to as the 

threshold region. The region above the last experimental measurement will be collectively 
denoted as the continuum region. The first experimental data close to the B meson threshold 
were taken by the CLEG [50-52] and CUSB [53] collaborations. The measurements at each 
c.m. energy have a 6% systematic uncertainty. More recently the BABAR collaboration [46] 
has measured the R^-ratio in the energy region between 10.54 GeV and 11.20 GeV, with 
significantly higher statistics and better control of systematic uncertainties (of the order of 
3%). These measurements are taken in small energy bins, densely populating the threshold 
region. The BABAR data supersedes the older data of CLEO and CUSB, and it has already 
been used in Refs. [27, 31, 41], in which the bottom mass was also determined. 

This BABAR data for the R^-ratio has not been corrected for initial-state radiation 
and vacuum polarization effects. Moreover, the effect of the T(45) resonance has not 
been subtracted,® so we have performed the subtraction ourselves, using the Breit-Wigner 
approximation and using for the total width the PDG value r 45 = 20.5 MeV: 


oBW/ N, ^ 9 M|g _ £45 _ 

a{MisY (s - -p TlgMjg ' 


( 6 . 2 ) 


For the subtraction of the T(45) resonance and the correction for the initial state 
radiation we take an approach similar to Ref. [27]. 


6.2.1 Subtraction of the T(4S') Radiative Tail 

Before subtracting the radiative tail of the T(45*) resonance one has to account for vacuum 
polarization effects. BaBar experimental data has been normalized to the theoretical Born 
level dimuon cross section (using the fine structure constant rather than the running effective 
electromagnetic coupling), instead of normalizing to the number of events with muons in 
the final state. Therefore one has to multiply the BaBar data with \a{s)/■, which we 
take as constant with value 0.93. 

The contribution to be subtracted from the BABAR data (already corrected for vacuum 
polarization effects) is the ISR-distorted tail of the T(45), which reaches to energies above 

®The radiative tails of the first three resonances are provided by BABAR, so they can be subtracted at 
the data level, before correcting for vacuum polarization effects. 
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its mass. The cross section R and the ISR-distorted cross-section R are related by a 
convolution relation 


— / — G{z,s)R{sz), (6-3) 

which can be used to determine the ISR effects on the T(4S') resonance given in 
Eq. (6.2). Here the lower integration bound is zq = (10 GeV)^/s. This value is not fully fixed 
by theoretical arguments, and it is chosen such that it excludes the narrow resonances, but 
keeps the major part of the T(4S') line shape. The radiator function G is given as [54, 55] 


G( 2 ;,s) = (1 - ^G{z,s), (6.4) 

G{z, s) = m eV.(s) F{s) [ 6^+^is) + Sg{s, z) ] , 

where the specific form of 13, F and the two (5’s can be found in Eq. (7) of [27]. Note 
that the function G{z, s) is divergent as z —>■ 1, but since 0 < /3 — 1 < —1, it is integrable. 
The divergent behavior is absent in G, which in the limit z —>■ 1 reduces to 

G(l,s) = /3(s) F(s) (5^+^(s). (6.5) 

After subtracting the radiative tail of the T(4S') we find that to a good approximation 
the cross section vanishes for energies below 10.62 GeV. Therefore we add an additional 
point to our BABAR dataset: i?fe(10.62 GeV) = 0 and take Rb = 0 for energies below 
10.62 GeV. Since the subtracted cross section does not exactly vanish between 10.5408 GeV 
and 10.62 GeV, we take the (small) contribution of the subtracted cross section in that 
region to the moments as an additional source of systematic correlated uncertainty. 

6.2.2 Deconvolution of Initial-State Radiation 

After subtraction of the radiative tails and correcting for vacuum polarization effects, the 
BABAR threshold data are corrected for ISR. The inversion of the convolution in Eq. (6.3), 
can be carried out in an iterative way [27]. Defining 5G{z, s) = G{z, s) — 5(1 — z) one can 
use a successive series of approximations 


/■l Ap’ 

R^{s) = R^{s)— / — 6G{z, s) R^~^{s z), (6-6) 

J Zo ^ 

where we denoted the j-th approximation of R{s) as R^ {s) and use as starting point 
R^{s) = R{s), the BABAR data after correcting for vacuum polarization effects and sub¬ 
tracting the radiative tails. In Eq. (6.6) we take zq = (10.62 GeV)^/s, using as a starting 
point the energy value for which the cross section vanishes after the subtraction of the ra¬ 
diative tails. To isolate the singularity at the higher endpoint one can perform a subtraction 
at z = 1, resulting in: 
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/Zo 


1 


/3(^) 


G(l,s)i?^-Hs)(l-^o) 


m 


(6.7) 


We use the trapezoidal rule to evaluate the integration on the discrete set of experi¬ 
mental data measurements labeled by the index i. Changing the integration variable from 
z to energy we find 


Ri=IiS + R{ 


( ’ ‘) ■ ( E?j ( E?-Ef 
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where we have used R\ = 72(10.62 GeV) = 0 for all iterations. After applying the 
procedure as many times as necessary to obtain a stable solution, one obtains the ISR- 
corrected cross section. Among the experimental measurements one finds two data points 
taken at very similar values of the energy: 10.86 GeV and 10.8605 GeV. It turns out that 
the fact that they lie very close makes the iterative procedure unstable. Therefore we drop 
the latter point from our analysis. 

In Fig. 8 we show the BABAR data after the subtraction of all radiative tails, before 
(red) and after (blue) ISR and vacuum polarization corrections. 


6.2.3 Determination of the Unfolding Error Matrix 

The BABAR collaboration splits the experimental uncertainties into statistical, systematic 
uncorrelated, and systematic correlated. We add the two former in quadrature to obtain 
the total uncorrelated uncertainty and rename the latter as the total correlated un¬ 

certainty The removal of the radiative tails of the T mesons has no effect on these 
uncertainties. Therefore, the correlation matrix for the BABAR data after the subtraction 
of the radiative tails, before it is corrected for ISR effects, can be written as 


MO/ = (e“)^ + er^r • (6.9) 

One needs to compute a new correlation matrix after each iteration. In this way we 
determine the unfolding error matrix. 

The master formula in Eq. (6.8) can be cast in a matrix form as follows: 


K = + (6-10) 

k=2 

where Rj is to be thought as the i-th component of the column vector R^, and Gik 
represents the (i, /c)-component of a matrix G. Here 72] depends only on the initial value 
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Figure 8. BABAR experimental data before (blue) and after (red) the ISR correction is applied. 
The purple bar on the right refers to the pQCD prediction for the continuum region. We have 
removed one data point at A = 10.8605 GeV. 


i?? and the result of the previous iteration The Gik do not depend on or the 

iteration step j. Therefore, for the error propagation one uses 


^ = 5ik , 
dRl 


dR^ 

^ 

— ^ik ? 


dRl 


( 6 . 11 ) 


both of them j-independent. We will denote with the correlation matrix among 
the entries of the vector R^ for a given iteration j. We also hnd it convenient to introduce 
the correlation matrix among R^ and R^, referred to as Finally we use the notation 

= (M-^ ®)^. We hnd for the correlation matrix after j iterations: 


mJJ = + G + G (6.12) 

where the elements of the matrix are given in Eq. (6.9). We hnd that after hve iter¬ 
ations the result has converged already to a level well below the experimental uncertainties. 
Our unfolded BaBar data agrees well with that worked out in Ref. [41]. 

6.2.4 Contribution of the Threshold Region 

After having corrected BABAR data for ISR and vacuum polarization effects, we use the 
trapezoidal rule for integrating the threshold region between 10.62 GeV and 11.20 GeV: 
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(6.13) 


Mf'' = — 

^ 2n 
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+ Rn 
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where Ri has been already ISR corrected and N is the number of data points. We have 
added the boundary condition point Ri = i?(10.62) = 0. From Eq. (6.13) one can compute 
the correlation matrix among for various n values, using the unfolding matrix among 
the Ri computed in Sec. 6.2.3. 


6.3 Continuum Region 

For the determination of the experimental moments from the region above 11.2 GeV we use 
pQCD (which has essentially negligible errors) supplemented by a modeling uncertainty. 
Comparing pQCD (purple line in Fig. 9) to a linear ht to the BaBar data in the region 
between 11.06 GeV and 11.2 GeV (red dotted line in Fig. 9) we hnd a 10% discrepancy con¬ 
cerning the central values. The ht function has a roughly constant 4% relative uncertainty. 
The ht linear function shows a growing pattern such that it would meet the pQCD predic¬ 
tion at around 11.5 GeV. This result is very robust, since a quadratic ht yields the same 
meeting point. To model the continuum in the region between 11.2 GeV and 11.52 GeV 
we patch together the linear ht function to the BaBar data and the result of pQCD above 
11.52 GeV using a cubic function, demanding continuity and smoothness at 11.2 GeV and 
11.52 GeV. The result is shown as the central red line in Fig. 9. Given that the relative 
discrepancy between experiment and pQCD for Rf, at the Z-pole is about 0.3% [56], we 
adopt a relative modeling error that decreases linearly from 4% at 11.2 GeV to 0.3% at mz, 
and stays constant for energies larger than mz- This is shown as the red band in Fig. 9. 
This uncertainty makes up for 96.9% of the total error for the hrst moment MY (which 
has an total 2.45% relative error), and 86.15% of the second moment MY (which has a 
total 1.85% relative error). Note that if we would adopt a constant 4% error for all energies 
above 11.2 GeV, this continuum uncertainty would make up for 97.24% of the total error 
for the first moment MY (from a total 2.60% relative error), and 86.46% of the second mo¬ 
ment MY (from a total 1.87% relative error). The difference is small because contributions 
from higher energies are suppressed. Following our procedure in Ref. [6] we consider this 
uncertainty as fully correlated for the various moments, but without any correlation to the 
narrow resonances or the threshold region. 

The perturbative QCD theoretical expression which we use to determine this contribu¬ 
tion includes the non-singlet massless quark cross section supplemented with bottom mass 
corrections up to 0(m^/s^).® It takes into account only contributions from the electromag¬ 
netic current coupled to the bottom quark. It reads: [57-61]^^ 

®We note that the double massive fermion bubble contribution to Rhb in Eq. (6.15) includes both virtual 
and real radiation terms in the large energy expansion. However, when this formula is used to compare 
pQCD to the existing BABAR data, below the four-bottom-quarks threshold, the real radiation should be 
excluded. We have checked that this inconsistency has an effect below 0.1%. 

^°The authors of Ref. [41] use the pole mass instead of the MS, and include and QED corrections. 
This explains some numerical differences in the analyses. 
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Figure 9. Comparison of ISR-corrected BABAR data in the continuum region (black dots with 
error bars) with pQCD (purple band). The red band shows our reconstruction of the continuum, 
which includes a linear fit to the BaBar data, patched to the pQCD prediction in a smooth way 
using a cubic polynomial in the energy. 


Qt m^{y/s),nf = 5 
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where 
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+ (^) (3545.81 - 158.311 + 9.66667- 538.813 + 37.8611 Ls 


-1037.79 L2 + 185.389 7.3) 


(6.15) 


with 


Ls = ln(^^^, Lm ^ In i ois = ■ ( 6 - 16 ) 

We use the initial conditions nihijnh) = 4.2 GeV and as{mz) = 0.118. 

Therefore, for the continuum region we use the following expression, 
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n 

Resonances 

10.62 - 11.2062 

11.2062 -00 

Total 

1 

1.394(12|22) 

0.270(2|9) 

2.862(0|108) 

4.526(12|111) 

2 

1.459(12|22) 

0.226(1|8) 

1.148(0|45) 

2.834(12|51) 

3 

1.538(12|22) 

0.190(1|7) 

0.611(0|24) 

2.338(12|34) 

4 

1.630(13|22) 

0.159(1|6) 

0.365(0|15) 

2.154(13|27) 


Table 3. Results for our computations of the experimental moments. The second column collects 
the contribution from the first four T resonances (using the narrow width approximation). The 
third to fifth columns show the contributions from the threshold (using ISR-corrected BABAR 
data) and continuum (using an interpolation between a linear fit to the BaBar data with highest 
energy and pQCD as a model for the lack of data) regions, and the total moment determinations, 
respectively. The two numbers quoted in parentheses correspond to the uncorrelated and correlated 
experimental uncertainties, respectively. All numbers are given in units of GeV”^''. 
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(6.17) 


with So = (11.2062 GeV)^, si = (11.52 GeV)^ and is a cubic function that 

smoothly interpolates between the linear fit to BaBar data and pQCD. Here 7 ' is the 
auxiliary variable used to parametrize our uncertainty, which we consider as 100 % correlated 
among the various moments. The related entries of the correlation matrix are trivially 
computed as 


PQCD ^ 

nn' Q^i 


^^pQCD 


(6.18) 


6.4 Final Results for the Experimental Moments 

The full result for the experimental moments is obtained by summing up all the portions 
described before, 


M“P = + MpQ^° . (6.19) 

We determine two correlation matrices among the first four moments. One of them 
comes from the various uncorrelated uncertainties, whereas the other encodes the systematic 
uncertainties. We denote them as the correlated and uncorrelated correlation matrices, 
respectively. These are computed by summing up the respective individual matrices from 
each region, and in the same way as we did for our charm analysis [ 6 ], we assume there is 
no region-to-region correlation. We find: 
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/oexp 

^uc 


/ 0.0002 0.0002 0.0002 0.0002 \ 
0.0002 0.0002 0.0002 0.0002 
0.0002 0.0002 0.0002 0.0002 
\ 0.0002 0.0002 0.0002 0.0002 J 


^exp 

^cor 


/ 0.0122 0.0055 0.0032 0.0021 \ 
0.0055 0.0026 0.0017 0.0012 
0.0032 0.0017 0.0011 0.0009 
\ 0.0021 0.0012 0.0009 0.0007 y 


where the (n, m) entry of each matrix is given in units of , 

and the total correlation matrix is the sum of C'uc’’ and Ctof ■ The contribution of each 
region to the final experimental moments and the corresponding uncertainties are presented 
in Table 3. 


7 Comparison to other Determinations of the Experimental Moments 

In this section we compare our result for the experimental moments for the bottom vector 
current correlator with previous determinations. These are collected in Table 4.^^ The most 
relevant comparison is between the second and third columns, where the most recent data 
on the narrow resonances and the BABAR continuum data are used. For the contributions 
from the narrow resonances we have a perfect agreement with [27], although slightly larger 
errors. For the threshold region our results are slightly smaller, and our uncertainties are 
almost identical; however, this is not a one-to-one comparison, since our integration region 
is slightly smaller. Indeed if we consider their energy range we agree with their numbers 
almost perfectly. The main difference between these two determinations is the estimate 
of the uncertainties coming from the continuum region, where the pQCD prediction for 
the Rft-ratio is employed. Whereas we adopt the more conservative approach described in 
Sec. 6.3, Ref. [27] employs only the perturbative uncertainties related to the purple band in 
Fig. 9. In Ref. [41] the same collaboration presents a more critical analysis of their errors. 
In particular they observe that the last experimental measurement of BABAR, after being 
corrected for ISR, disagrees with the pQCD prediction at the 20% level (way outside the cor¬ 
responding uncertainties).^^ To resolve this discrepancy they assume two possible scenarios: 
a) pQCD starts being reliable at energies above 13GeV (therefore the authors interpolate 
between the last experimental point and pQCD at 13GeV); b) BABAR systematic errors 
have been underestimated (therefore the central values of the experimental measurements 
are rescaled by a factor of 1.21). Ref. [41] quotes the values of the experimental moments 
and the resulting values for the bottom mass for these two scenarios. Since the effect of 
these differences of the two bottom masses obtained from is only slightly larger than 
the size of the other uncertainties (that is, the uncertainties of the theoretical moments 
plus the other experimental errors) added quadratically, it is argued that this issue can 

the case of Ref. [7], we reconstruct the experimental moments from their Table 3, where the moments 
are split in several different contributions. For the reconstructed uncertainty, we take one half of the error 
of the narrow resonances correlated to each other, and the other half as uncorrelated. The errors from 
patches where theory input is used are taken as fully correlated to one another. The total narrow-resonance 
error, and the total “theory-patch” error are added in quadrature to get the final uncertainty. 

^^From our own computation of the ISR-corrected 7?i,-ratio, we only observe a 10% deviation between 
the last data point and the pQCD prediction, see Fig. 9. 
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be ignored. We disagree with this argument, since the issue constitutes an independent 
source of uncertainty not covered by the other errors and, in particular, being unrelated to 
uncertainties in the theoretical moments. Therefore this shift must be taken as an addi¬ 
tional source of error on the experimental moments (and indeed would then dominate the 
corresponding total error). The additional error (to be added in quadrature to the one in 
round brackets) is quoted in square brackets in the third column of Table 4. It amounts 
to an additional error of 30,18,11 and 7MeV for fnh{fni,) extracted from moments MY to 
MY , respectively. 

Refs. [7, 17, 44] have used the older CLEO and CUSB experimental measurements, 
resulting in relatively large uncertainties. In Ref. [44] the CLEO measurements are divided 
by a factor of 1.28, and an error of 10% is assigned. It is argued that this procedure is 
necessary to reconcile old and new CLEO measurements, as well as to improve the agreement 
with pQCD predictions. Ref. [7] uses values for the T-states electronic partial widths given 
by the PDG 2002, which have larger uncertainties. This makes their determination of the 
experimental moments rather imprecise. 

Concerning the continuum region where no measurements exist, while some previous 
analyses have taken a less conservative approach than ours, in Ref. [62] a much more 
conservative approach is adopted. In this region they consider the i?;,-ratio as constant 
with a 66% uncertainty. In Ref. [7] also a more conservative approach is adopted. Between 
11.1 and 12GeV 0{a'^) pQCD errors are used, which are larger than 10%; for energies 
above 12 GeV a global 10% correlated error is assigned. 


n 

This work 

Chetyrkin et al. ’09 [27] 

Kuhn et al. ’07 [44] 

Corcella et al. ’03 [7] 

1 

4.526(12|111) 

4.592(31)[67] 

4.601(43) 

4.46(17) 

2 

2.834(12|51) 

2.872(28)[51] 

2.881(37) 

2.76(15) 

3 

2.338(12|34) 

2.362(26) [40] 

2.370(34) 

2.26(13) 

4 

2.154(13|27) 

2.170(26)[35] 

2.178(32) 

2.08(12) 


Table 4. Comparison of our results for the experimental moments of the bottom vector correlator 
(2nd column) to previous determinations (3rd to 5th columns). The 2nd and 3rd columns use 
BABAR data from Ref. [46], while 4th and 5th use older data from Refs. [50, 51[. The 3rd and 
4th columns use perturbative uncertainties in the continuum region, while 2nd and 5th use a more 
conservative estimate based on the agreement of data and pQCD. In the 3rd column, we quote in 
square brackets our own estimate of an additional systematic error from the considerations made 
in Ref. [41[. All numbers are given in units of GeV~^”. 


8 Computation of the Experimental Values for the Ratios of Moments 

Once the experimental values for the moments of the vector and pseudoscalar correlators 
have been computed, it is in principle a straightforward exercise to calculate ratios of 
them. The central value is obtained by simply taking the ratio of the corresponding central 
values. To obtain the uncertainties (or more generally, the correlation matrix among the 
different ratios of moments), one needs to have access to the complete correlation matrix 
among the moments. Our computation in Ref. [6] [see Eqs. (3.21) and (3.22)] for the 
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Figure 10. Comparison of various determinations of the experimental moments for the bottom 
vector correlator. Results in blue correspond to analyses of the same collaboration. The green 
result and the determination at the top do not use the new BABAR results. 


charm experimental moments, and the procedure presented in Sec. 6 [see Eq. (6.20)] to 
determine the bottom experimental moments, yield the two desired correlation matrices, for 
statistical and systematical correlations. For the pseudoscalar moments the information on 
correlations is not provided in Ref. [8]. Therefore we make the simplest possible assumption, 
which is that the moments are fully uncorrelated. This will most certainly overestimate 
the uncertainties for the ratios of moments, but given that we are anyway dominated by 
perturbative uncertainties, our approach appears justified. We collect our results for the 
computation of the ratios of experimental moments in Table. 5. Readers interested in the 
full correlation matrix among them can send a request to the authors. 

9 Results 

In this section we present the final results for our analyses at O(a^). We take method (c) 
(linearized iterative expansion) as our default expansion. For the estimate of the perturba¬ 
tive uncertainty, we perform double scale variation in the ranges rficijnc) < < 4 GeV 

for charm (either correlator), and mb{fnb) < fXa, < 15 GeV for bottom, and we discard 
3% of the series with the worst convergence (that is, with highest values of the 14 con¬ 
vergence parameter). For the charm mass determinations (either vector or pseudoscalar 
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n 

Vector Uf = 4: 

Vector Uf = 5 

Pseudoscalar Uf = 4 

1 

6.969(32|59) 

6.262(10|53) 

0.971(32) 

2 

8.807(23|26) 

8.251(09|48) 

1.048(53) 

3 

9.547(14|13) 

9.212(08|35) 

1.092(77) 


Table 5. Ratios of experimental moments for the vector correlator with 4 and 5 flavors (second and 
third column, respectively), and for the pseudoscalar correlator with 4 flavors (fourth column). For 
the vector current, the first error in parenthesis corresponds to the statistical uncertainty, whereas 
the second corresponds to the systematic one. For the pseudoscalar correlator we only quote the 
lattice error. Moments given in units of 10“^ GeV“^, 10“^ GeV“^ and, 10“^ GeV“^ for the second, 
third, and fourth column, respectively. 


correlator) we use the first moment as our default, given that it is theoretically more reli¬ 
able than the higher moments. For the analysis of the bottom mass from regular moments, 
we use as our default, since it is less afflicted by systematic experimental errors than 
the first moment, and is nevertheless theoretically sound. For the charm and the bottom 
mass analyses we also examine the ratio of the second over the first moment as a cross 
check and validation of the results from regular moments. The results for the experimental 
moments are collected in: the last column of Table 9 in Ref. [6] (charm vector correlator 
regular moments); the last column of Table 3 (bottom vector correlator regular moments); 
Table 1 (lattice regular moments); and Table 5 (all ratios of moments). 

We also analyzed higher (and also lower for the case of bottom) moments and ratios of 
moments for all correlator and quark species. Since, as already discussed, fixed-order and 
contour-improved higher moments are particularly afflicted by their nonlinear dependence 
on the quark mass, we only consider the linearized and iterative methods for this analysis. In 
any case, since higher moments have a larger sensitivity to infrared effects and are therefore 
theoretically less sound, the analysis involving higher moments mainly aims at providing 
cross checks. The results are collected in a graphical form in Fig. 11, and the numerical 
results can be obtained from the authors upon request. 

Our final determinations include nonperturbative effects through the gluon condensate 
including its Wilson coefficients at order 0{as). Furthermore, we assign as a conservative 
estimate of the nonperturbative uncertainty twice the shift caused by including the gluon 
condensate. In any case, this error is very small, particularly for the bottom mass determi¬ 
nation. One source of uncertainty which we have not discussed so far is that coming from 
the strong coupling constant. Although the world average as{mz) = 0.1185 ± 0.006 has 
a very small error, see Ref. [48], one cannot ignore the fact that it is fairly dominated by 
lattice determinations, e.g. [29]. Furthermore, there are other precise determinations with 
lower central values and in disagreement with the world average from event-shapes [63-66] 
and DIS [67]. A review on recent Ug determinations can be found in Refs. [68-70]. There¬ 
fore, in analogy with Ref. [6], we perform our analyses for several values of as{mz) between 
0.113 and 0.119, and provide the central values and perturbative errors as (approximate) 
linear functions of as{mz)- The other uncertainties are essentially cis-independent, so we 
just provide the average. We also quote quark mass results for Og taken from the world 


- 31 






n n n 


(a) 


(b) 


(c) 



n n n 


(d) 


(e) 


(f) 


Figure 11. Charm and bottom quark mass determinations for different moments (upper row) or 
ratios of moments (lower row), for the linearized (in blue) and iterative (in red) methods. Panels 
(a), (b) [(e), (f)] show the results for the charm mass from moments [ratios of moments[ of the 
vector and pseudoscalar correlator, respectively. Panels (c) and (g) show the results for the bottom 
mass from the vector correlator, for moments and ratios of moments, respectively. 


average: 


a^(mz) = 0.1185 ±0.0021, (9.1) 

where we adopt an uncertainty 3.5 times larger than the current world average [48]. 
We note that in Ref. [6] we have taken as{mz) = 0.1184± 0.0021 as an input which causes 
only tiny sub-MeV differences in the quark masses. We refrain ourselves from presenting 
the as dependence of the higher-moment result, which the reader can get from the authors 
upon request. 

For the numerical analyses that we carry out in this article we have created two com¬ 
pletely independent codes: one using Mathematica [71] and another using Fortran [72], 
which is much faster and suitable for parallelized runs on computer clusters. The two codes 
agree for the extracted quark masses at the level of 1 eV. 

9.1 Results for the Charm Mass from the Vector Correlator 

For the analysis using the first moment of the charm vector correlator we use the experimen¬ 
tal value quoted in Eq. (4.1) of Ref. [6]: = (0.2121 ± 0.0020stat ± O.OOSOgyst) GeV~^. 

The outcome of this analysis, and one of the main results of this paper, is: 
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Figure 12. Dependence on as{mz) of the central values of md^c) o-i^d the corresponding per¬ 
turbative (red), statistical (orange), systematical (blue) and nonperturbative uncertainties (green), 
for the analysis of the first moment [panels (a) and (b)] and the ratio of the second over the first 
moment, [panels (c) and (d)[, corresponding to the vector correlator. 


mdruc) =1.288 ± (0.006)stat ± (0.009)syst ± (0.014)pert (9.2) 

± (0.010)„, ± (0.002) (gg> GeV, 

where the quoted errors are (from left to right) experimental uncorrelated, experimental 
correlated, peturbative, due to the uncertainty in as as given in Eq. (9.1), and nonpertur¬ 
bative. If we adopt the correlated scale variation 2 GeV < f^a = < 4 GeV, we obtain for 

method (c) 1.297 ± (0.005)pert) with the other errors essentially unchanged. For method (a) 
we would get 1.290 ± (0.0007)pert) with a scale variation even smaller than the nonpertur¬ 
bative uncertainty, and 20 times smaller than our perturbative error estimate with double 
scale variation [3 times for method (c)]. The dependence on as{mz) is shown graphically 
in Figs. 12(a) and 12(b), and analytically the result reads: 

mdmc) = (1.288 + 4.40 x [a.(mz) - 0.1185]) ± (0.006)stat ± (0.009)syst (9.3) 

± (0.014 + 0.95 X [as{mz) - 0.1185])pert ± (0.002 )(gg> • 
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Eqs. (9.2) and (9.3) supersede the results given in Eqs. (4.5) and (4.2) of Ref. [6], 
respectively. 

For the ratio of the second over the first moment of the vector correlator we use as the 
experimental input = (6.969 ± 0.032stat i 0.059syst) x 10“^GeV“^, which yields 

the following result for the charm mass: 

mc(mc) =1.271 ± (0.003)stat ± (0.005)syst ± (0.016)pert (9.4) 

± (0.004)„^ ± (0.004) (GG>GeV. 

With correlated variation 2 GeV < //« = Hm < 4GeV we get 1.258 ± (0.005)pert and 
1.279 ± (0.007)pert for methods (a) and (c), respectively. In this case the scale variations 
are a factor of 2 to 3 smaller than our perturbative error estimate. The q;<j dependence, 
which can be seen in Figs. 12(c) and 12(d), has the form: 


mc(m,) = (1.271 + 1.64 x K(mz) - 0.1185]) ± (0.003)stat ± (0.005)syst (9.5) 

± (0.016 + 1.081 X [as{mz) - 0.1185])pert ± (0.004 )(gg> • 

We observe that the central value for the ratios of moments is 17 MeV smaller than 
for the first moment analysis, but fully compatible within theoretical uncertainties. Fur¬ 
thermore, the dependence on as of the central value obtained from the regular moment 
analysis is larger, which translates into a corresponding larger error due to the uncertainty 
in as- Both determinations have very similar perturbative uncertainties for any value of 
as- We also see that the charm mass from the ratio of moments has smaller experimental 
uncertainties, as a result of cancellations between correlated errors. Moreover, the charm 
mass result from has a nonperturbative error twice as large as that from MY - The two 
charm mass results from the first moment and the moment ratio are compared graphically 
in Fig. 15(a). 

9.2 Results for the Charm Mass from the Pseudoscalar Correlator 

For the analysis of the first moment of the charm pseudoscalar correlator we employ 
j^P,\att _ ^g_]^ 4 g 2 0.0020iatt) GeV“^ [8], which yields the following charm mass de¬ 

termination: 


mc(mc) = 1.267 ± (0.008)iat ± (0.035)pert ± (0.019)„, ± (0.002 )(gg> GeV. (9.6) 

With correlated scale variation 2 GeV < fJ-a = fJ-m < 4 GeV we obtain the central values 
1.278 and 1.276 GeV, for methods (b) and (c), respectively. In both cases the scale variation 

^®Had we taken the fixed-order expansion (a) and correlated scale variation 2 GeV < fJ-a = Pm < 4 GeV 
as the estimate for the perturbative uncertainty, the result from BX with all errors added quadratically 
would be 1.258 ± 0.013 GeV, whereas the result from MY would read 1.290 ± 0.015 GeV. Both results would 
not be consistent to each other. 
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Figure 13. Dependence on asijnz) of the central values of mc(jnc) and the corresponding 
perturbative (red), lattice (blue) and nonperturbative uncertainties (green), for the analysis of the 
first moment [panels (a) and (b)] and the ratio of the second over the first moment, [panels (c) and 
(d)[, corresponding to the pseudoscalar correlator. 


is 4 MeV, 8 times smaller than our perturbative error estimate with double scale variation. 
For the ag dependence, we hnd 


rricimc) = (1.267 + 8.36 x [a^(mz) - 0.1185]) ± (0.008)iat (9.7) 

± (0.035 + 2.38 X [as{mz) - 0.1185])pert ± (0.002 )(gg} , 

which is also displayed in Figs. 13(a) and 13(b). As expected, the perturbative error is 
much larger than for the vector correlator, and has a stronger dependence on ag. We see 
that the central value has a much stronger dependence on ag as well, which again translates 
into a large error due to the uncertainty in the strong coupling. The central value is 21 MeV 
lower than Eq. (9.2), but fully compatible within errors [see Fig 15(a)]. The nonperturbative 
uncertainties are identical to the vector current case. 

For the ratio of second over the hrst moment of the pseudoscalar correlator we use 
^P, latt _ 0.0032iatt) GeV“^. We hnd for the charm mass 
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mc(me) = 1.266 ± (0.020)iatt ± (0.018)pert ± (0.006)„, ± (0.002 )(gg> GeV. (9.8) 


Using correlated variation 2GeV < < 4GeV one obtains 1.270 ± (0.007)pert 

and 1.278 ± (0.003)pert for methods (a) and (c), respectively. These scale variations are a 
factor 3 and 6 smaller than our perturbative error estimate, respectively. The Ug dependence 
is 


mc(mc) = (1.266 + 2.31 x [a^(mz) - 0.1185]) ± (0.020)iatt (9.9) 

± (0.018 + 1.25 X [as{mz) - 0.1185])pert ± (0.002 )(gg> • 

The central values for both Mf and Rf are almost identical, but their ag dependence 
is not: the latter is much smaller (even smaller than for MY, but larger than for RY)- Note 
that the lattice error is larger for the ratio since we made the very conservative assumption 
that they are fully uncorrelated. This is because the correlation matrix for various lattice 
moments is unknown. The perturbative error reduces by a factor of two for any value 
of ag when using the ratio, but we have checked that this only happens for the iterative 
expansion. On the other hand, the Ug dependence of the perturbative uncertainty is smaller 
for the regular moment determination. The nonperturbative errors are identical. 

All charm determinations are illustrated graphically in Fig. 15(a), where in red we show 
our preferred determination from the first moment of the vector correlator. 

9.3 Results for the Bottom Mass from the Vector Correlator 

For our determination of the bottom quark mass from the second moment of the vector 
correlator we use for the experimental moment = (2.834 ± 0.012stat =t 0.051syst) x 

10“^ GeV“^, and we obtain 


mfc(m6) =4.176 ± (0.004)stat ± (0.019)syst ± (O.OlO)pert (9.10) 

± (0.007)„, ± (0.0001)(GG>GeV. 

The perturbative error is 30% smaller than for the charm vector correlator analysis, as 
a result of the smaller value of q;<j at the scales close to the bottom mass. This is consistent 
with our discussion on the convergence properties of perturbation series for the bottom 
quark carried out in Sec. 4. The total error is dominated by the experimental systematic 
uncertainty, which in turn mainly comes from the continuum region where one relies on 
modeling in the absence of any experimental measurements. The nonperturbative error is 
completely negligible. This is expected since it is suppressed by two powers of the bottom 
mass. Using the correlated scale variation 5 GeV < Ha = < 15 GeV for methods (a) 

and (c) we get 4.178 and 4.182 for the central values and scale variations which are 20 and 
3 times smaller, respectively. The ag dependence reads 
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Figure 14. Dependence on as(™z) of the central values of mi,{Wii,) and the corresponding per¬ 
turbative (red), statistical (orange), systematical (blue) and nonperturbative uncertainties (green), 
for for the analysis of the second moment [panels (a) and (b)] and the ratio of the second over the 
the first moment [panels (c) and (d)[ of the vector correlator. 


fnbimb) = (4.176 + 3.22 x [a,(mz) - 0.1185]) ± (0.004)stat ± (0.019)syst (9.11) 

± (0.010 + 0.472 X [asimz) - 0.1185])pert ± ( 0 . 0001 )(gg> • 

For the analysis based on the ratio of the second moment over the first, we use = 

(6.262 ± O.OlOstat i 0.053syst) x 10 “^GeV“^, and with this value we obtain the bottom 
mass 


mb{mb) =4.179 ± (0.003)stat ± (0.017),yst ± (0.009)pert (9.12) 

± (0.003)„, ± (0.0002) (gg> GeV. 

With correlated scale variation 5 GeV < fJ-a = fJ-m < 15 GeV we obtain 4.175 ± 
(0.003)pert and 4.182 ± (0.004)pert for methods (a) and (c), respectively. In this case 
the scale variation is smaller by a factor 3 and 2, respectively. The as dependence reads as 
follows: 
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Figure 15. Charm (a) [bottom (b)] mass determinations from the first [second] moment of the 
vector correlator (in red), the first moment of the pseudoscalar correlator (green, charm only), 
and the ratio of the second over the first moment of the vector (blue) and pseudoscalar correlator 
(purple, charm only). 
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mb(mfe) = (4.179 + 1.199 x [asimz) - 0.1185]) ± (0.003)stat ± (0.017)syst (9.13) 

± (0.009 + 0.426 X [as(mz) -0.1185])pert ± (0.0002 )(gg> • 

Although the central value for the ratio analysis is 3 MeV higher, this has no significance 
given the size of the uncertainties. The dependence of the central value on as is three times 
smaller for the ratio analysis. The perturbative error and its as dependence are roughly 
the same for the ratio and the single moment analysis. Moreover, the two experimental 
errors are very similar. This is because, even though there is some cancellation of correlated 
errors in the ratio, a signihcant part of the huge systematic error of the first moment remains 
uncanceled. 

A graphical illustration of the two bottom mass determinations is shown in Fig. 15(b). 
Both combined uncertainties and central values are rather similar, and we adopt the result 
from the second moment (in red) as our default result. 

10 Comparison to other Determinations 

In this section we make a comparison to previous analyses of our updated charm mass 
determination from the vector correlator, our new results for the charm mass from the 
pseudoscalar current correlator and of our bottom mass determination. We restrict our 
discussion to determinations which use QCD sum rules with infinite as well as hnite energy 
range for the vector or pseudoscalar current correlators, and including relativistic and non- 
relativistic versions of the sum rules. We do not cover charm mass determinations from DIS 
or bottom mass determinations from jets (which are in any case rather imprecise), as well 
as determinations which are based on the mass of bound states (B mesons or quarkonia) or 
B decays. 
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In Figs. 16(a) and 16(b) we present in a graphical form a compilation of recent sum 
rule determinations of the charm and bottom masses, respectively. We have labeled them 
from top to bottom with numbers from 1 to 14. We note that comparing these results, 
one has to keep in mind that different analyses in general employed different values and 
uncertainties for the strong coupling. Only the analyses in Refs. [6, 27, 38, 47, 73] and ours 
have provided the dependence of their results on the value of as{mz)- 

10.1 Charm Mass 

Let us first focus our attention on the charm mass. Fig. 16(a). Within each color, de¬ 
terminations are ordered according to publication date. In red (determinations 12 to 14) 
we show the results of our collaboration: 12 and 13 for the vector correlator, the former 
(dashed) corresponding to Ref. [6] without trimming procedure, and the latter (solid) cor¬ 
responding to this work, which includes the trimming procedure. Determinations 12 to 14 
are the only analyses using uncorrelated scale variation. Determination 6 (gray) [38] sets 
^rn = T^cijnc)) and all the other analyses have set = /^a- Determinations in blue (1 [29], 
2 [28] and 3 [8]) were performed by the HPQCD collaboration, which employ method (b) 
for the pseudoscalar correlator moments used for the mass determination in their lattice 
analyses. Only 1 to 3 and 14 use pseudoscalar moments, while all the other analyses use the 
vector correlator. Among those 7 [26], 8 [27], 9 [44] and 10 [17] use data in the threshold 
region only up to 4.8 GeV; analysis 6 uses two patches of data in the threshold region, one 
from threshold to 4.7 GeV, and another between 7.2 and 11 GeV; analyses 12 and 13 use 
all available data (see Ref. [6] for the complete bibliographic information on charm data). 
The result in 6 uses only 0{a‘l) perturbative input [all the other analyses utilize O^a^) 
computations] and older information on the narrow resonances, from the PDG 2006. They 
also study the hxed-order expansion and two methods of contour improvement. In black (8 
to 10) we display results using fixed order analyses at ©(af) from the Karlsruhe (8 and 9) 
and Wiirzburg (10) collaborations. 

Analysis 7 (orange) corresponds to weighted hnite-energy sum rules. They employ a 
kernel which enhances the sensitivity to the charm mass, and at the same time reduces the 
sensitivity to the continuum region. Green color analyses, collected in 4 [74, 75] and 5 [73], 
apply other kinds of sum rules. Analysis 5 uses a hnite energy sum rule similar to 7, but 
the kernel makes the sensitivity to the charm mass quite small. On the other hand, the 
two determinations of 4 use shifted moments, ratios of shifted moments, and exponential 
sum rules, and consider only the contributions from the first 6 vector resonances in the 
narrow width approximation, and the pQGD prediction for the continuum. The lower 
analysis of 4 [74] includes contributions from condensates up to dimension 6; the higher [75] 
includes condensates up to dimension 8. In purple (11 [76]) we show the only analysis 
which uses large-n moments for the charm mass hts, employing NRQGD methods to sum 
up large logs and threshold enhanced perturbative corrections, supplemented with hxed 
order predictions for the formally power suppressed terms. This analysis uses contributions 
from narrow resonances, plus a crude model for the threshold and continuum patches, for 
which a conservative uncertainty is assigned. We note, however, that this analysis might 


39 - 


be questioned since perturbative NRQCD is in general not applicable for the charmonium 
states. 

Our new vector correlator result agrees well with the world average, having a similar un¬ 
certainty. Our result is fully compatible with the other determinations shown in Fig. 16(a). 
As mentioned already before, we disagree with the small perturbative uncertainties related 
to the scale variations of the vector and/or pseudoscalar moments adopted in analyses 1 to 
3 and 7 to 10. 

10.2 Bottom Mass 

Let us now turn our attention to the bottom mass results, see Fig. 16(b). The coloring and 
chronological conventions are analogous to Fig. 16(a), and we try to keep a similar ordering. 
We show three nonrelativistic determinations (11 [47], 12 [62] and 13 [77]) in purple; 0(ag) 
fixed-order analyses are shown in gray (5 [7]), black (10 [78]) and green (4 [79]); finite energy 
sum rules also based on fixed-order appear in orange (6 [31]) and green (4); there are two 
lattice analyses in blue, collected in 1 [28, 30]. Analyses 3, 6, 7 and 11 to 13 use the new 
BABAR data, whereas the others use the older CLEO and CUSB data. Analyses 4 and 
11 include only the contributions of the first six vector resonances. Analyses 4 and 5 use 
older measurements of the electronic width for the narrow resonances. Analyses 3, 4, 6 to 
9 use pQCD in the high-energy spectrum for the experimental moments. The theoretical 
treatment of the bottom mass analyses in red, gray, black, blue and green are in complete 
analogy to their charm mass analyses: 3 and 4 for bottom correspond to 4 and 5 for charm, 
respectively; 1 in bottom corresponds to 1 and 2 for charm. 

The upper analysis of 1 [30] uses a nonrelativistic lattice action to compute ratios of 
large-n moments, which are later compared to relativistic continuum perturbation theory. 
Because the continuum computation do not sum up Sommerfeld enhanced terms, this proce¬ 
dure is questionable. The analysis 10 uses a combination of Mq with the infinite momentum 
transfer moment, both in fixed-order, in order to constrain the continuum region. They only 
use experimental information on narrow resonances, and model the rest of the spectrum 
with theory predictions. Finally, they make the following scale choice: [ia = = thb(mb) 

and estimate the truncation error from an ansatz for the 0(af) term.^'^ Analyses 11 and 13 
use large-n moments and NRQCD methods for their theoretical moments. Analysis 12 uses 
NRQCD fixed-order perturbation theory at N^LO (which accounts for the summation of the 
Coulomb singularities) and 13 uses renormalization group improved perturbation theory in 
the framework of vNRQCD at N^LL order. Both analyses employ low-scale short-distance 
masses to avoid ambiguities related to the pole mass renormalon. Analysis 12 also uses 
N^LO NRQCD fixed-order input, but is incomplete concerning the contributions from the 
continuum region in the theoretical moments. Moreover, they extract the pole mass which 
is then converted to the MS scheme. 

Our result 14 is in full agreement with the world average, having a slightly smaller 
uncertainty. It also agrees with the other analyses shown, with slightly smaller or compa- 

^^Ref. [78] also makes a determination of the charm mass. We exclude it from our comparison since it is 
not used in the PDG average. 
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rable uncertainties. We disagree with the small perturbative uncertainties related to scale 
variations quoted in 6 to 10. 
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Figure 16. Comparison of recent determinations of charm (a) and bottom (b) quark masses from 
sum rule analyses. Red results correspond to our determination. Black and gray correspond to 
©(ofg) and 0{ag) analyses, respectively. Purple results use nonrelativistic sum rules. Orange use 
weighted finite energy sum rules. Blue results are based on QCD sum rules using lattice simulation 
results as experimental data. Green labels other kinds of sum rule analyses (FESR, Q^-dependent 
moments, ratios of moments). 


11 Conclusions 

In this work we have determined the MS charm and bottom quark masses from quarko- 
nium sum rules in the framework of the OPE, using perturbative computations, 

plus nonperturbative effects from the gluon condensate including its Wilson coefficient at 
0(as). For the determination of the perturbative uncertainties we independently varied 
the renormalization scales of the strong coupling and the quark masses, in order to account 
for the variations due to different possible types of a* expansions, as suggested earlier in 
Ref. [6]. 

In order to avoid a possible overestimate of the perturbative uncertainties, coming 
from the double scale variation in connection with a low scale of as and resulting in badly 
convergent series, we have re-examined the charm mass determination from charmonium 
sum rules (vector correlator) supplementing the analysis with a convergence test. The 
convergence test is based on Cauchy’s radical test, which is adapted to the situation in 
which only a few terms of the series are known, and quantihes the convergence rate of each 
series by the parameter 14. We hnd that the distribution of the convergence parameter 14 
coming from the complete set of series peaks around its mean value, and allows to quantify 
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the overall convergence rate of the set of series for each moment in a meaningful way. This 
justifies discarding (or “trimming”) series with values of 14 much larger than the average. 
For our analysis we discard 3% of the series having the highest 14 values, which results in 
a reduction of the perturbative uncertainty for the MS charm mass rndrnc) from 19MeV 
in [6] to 14MeV (which amounts to 26%), and a small shift of + 5MeV in the central value. 
Our new determination of the charm mass from the first moment (which is theoretically 
the cleanest) of the vector correlator reads: 

Tficifnc) = 1.288 ± 0.020 GeV, [Vector Correlator] (H-l) 

where all sources of uncertainty have been added in quadrature. This result supersedes 
our corresponding earlier result from Ref. [6], which was 1.282 ± 0.024 GeV. This makes it 
clear that the trimming procedure discards series which produce small values of the charm 
mass. 

We have applied the same method of theory uncertainty estimate to analyze the 
HPQCD lattice simulation results for the pseudoscalar correlator. Our convergence test 
signals that the pseudoscalar moments have far worse convergence than the corresponding 
vector ones. This translates into an uncertainties of 35 MeV due to the truncation of the 
perturbative series and the error in Ug (roughly twice as big as for the vector determina¬ 
tion). In contrast, using correlated scale variation (e.g. setting the scales in the mass and 
the strong coupling equal) the scale variation can be smaller by a factor of 8. Our new 
determination from the first moment (again being the most reliable theoretical prediction) 
of the pseudoscalar correlator reads mc(mc) = 1.267 ±0.041 GeV, where again all individual 
errors have been added in quadrature. The combined total error is twice as big as for the 
vector correlator, and therefore we consider it as a validation of Eq. (11.1) in connection 
with the convergence test. The result is in sharp contrast with the analyses carried out 
by the HPQCD collaboration [8, 28, 29], where perturbative uncertainties of 4 MeV are 
claimed. We have checked that, as for the vector correlator, for the different possible types 
of CKs-expansion the correlated variation in general leads to a bad order-by-order convergence 
of the charm mass determination. 

The second important result of this work is the determination of the bottom quark 
mass from the vector correlator. We have reanalyzed the experimental moments by com¬ 
bining experimental measurements of the first four narrow resonances, the threshold region 
covered by BABAR, and a theoretical model for the continuum. This theoretical model 
is an interpolation between a linear fit to the BaBar data points with highest energy and 
pQCD, to which we assign a 4% systematic uncertainty which decreases linearly to reach 
0.3% at the Z-pole, and stays constant at 0.3% for higher energies. Our treatment is mo¬ 
tivated by the error function yielded by the fit to BaBar data in the energy range between 
11.0 and 11.2 GeV and the discrepancy between pQGD and experimental measurements at 
the Z-pole. This results into a large error for the first moment, and therefore we choose the 
second moment (which is theoretically as clean as the first one for the case of the bottom 
quark) for our final analysis, giving a total experimental uncertainty of 18 MeV. Our treat¬ 
ment of the experimental continuum uncertainty is in contrast to Ref. [27], where instead 
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the very small perturbative QCD uncertainties (less than 1%) are used, claiming an experi¬ 
mental uncertainties of 6 MeV. In the light of the analysis carried out here, supported by the 
observations made in Ref. [41], we believe this is not justified. Our convergence test reveals 
that, as expected for the heavier bottom quark, the perturbative series converge faster than 
for the charm quark. Correspondingly, the perturbative and Ug uncertainties are ~ 30% 
smaller than those for charm. Taking correlated scale variation as used in Refs. [27, 31] the 
perturbative error estimate can shrink up to a factor of 20. We also find that correlated 
variation leads to incompatible results for the different types of a^-expansions. Our final 
result for the bottom mass from the second moment, with all errors added in quadrature, 
reads: 


mbimb) = 4.176 ± 0.023 GeV, [Vector Correlator] (11-2) 

where the total error is fairly dominated by the systematic error, which comes from 
the continuum region of the spectrum. Our uncertainty is very similar to the one obtained 
by the HPQCD analysis, but 30% larger than the 16 MeV claimed by [27]. Our central 
value is 13 MeV larger than the latter. This good agreement is a result of two effects that 
push in opposite directions: smaller value of the second experimental moment, and different 
perturbative analysis. Curiously enough, a similar accidental cancellation was observed for 
the charm mass in [6]. 

In order to further validate the results discussed above, we have also analyzed the ratios 
of consecutive moments of each one of the three correlators as alternative observables. In 
all cases the results from the moment ratios agree very well the regular moment analyses. 
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A Numerical Values for the Perturbative Coefficients 

In this appendix we succinctly collect the numerical values for all of the coefficients appear¬ 
ing in the perturbative series, such that our analysis can be reproduced. We organize these 
values in tables, each of them corresponding to a different equation. 


43 - 



= 1 

n = 2 

n = 3 

n = 4 

[avirif 

O 

"S' 

II 

-4.011 

- 6.6842 

-9.7224 

- 13.0879 

[avirif 

o 

"S' 

II 

- 36.9604 

-69.9123 

-112.669 

- 165.326 

[ap{nf 

II 

o 

8.02101 

0 

-9.72244 

- 20.9406 

[ap{nf 

o 

II 

39.1439 

- 36.3842 

-152.67 

-309.925 


Table 6 . Numerical values for the coefficients of Eq. (2.10) for the vector correlator with n/ = 5 
(first two columns), and for the pseudoscalar correlator with nj = 4 (last two columns). (Gluon 
condensate contribution). 
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Table 7. Numerical values of the coefficients for Eq. (2.8) for 11^(0) in the MS scheme and nj = 5, 
and P{q^ = 0) for 71/ = 4. 










[Cvt^ 


n = \ 

i = 0 

1.06667 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

2.55473 

2.13333 

0 

0 

0 

0 

0 

0 

0 

i = 2 

3.15899 

8.33909 

4.17778 

0 

-4.89657 

-4.08889 

0 

0 

0 

i = 3 

- 7.76244 

18.2235 

29.3221 

8.12346 

- 18.2834 

-37.1221 

-16.0148 

9.38509 

7.83704 

n = 2 

i = 0 

0.457143 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

1.10956 

1.82857 

0 

0 

0 

0 

0 

0 

0 

i = 2 

3.23193 

7.20649 

5.40952 

0 

-2.12665 

- 3.50476 

0 

0 

0 

i = 3 

- 2.64381 

19.0805 

35.2229 

14.1249 

- 15.0705 

-32.0439 

- 20.7365 

4.07609 

6.71746 

n = 3 

i = 0 

0.270899 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.519396 

1.6254 

0 

0 

0 

0 

0 

0 

0 

i = 2 

2.06768 

5.57705 

6.43386 

0 

-0.995509 

-3.11534 

0 

0 

0 

i = 3 

- 1.17449 

14.8309 

36.8953 

21.0888 

-9.18132 

-25.3067 

-24.6631 

1.90806 

5.97108 

n = 4 

i = 0 

0.184704 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.203121 

1.47763 

0 

0 

0 

0 

0 

0 

0 

i = 2 

1.22039 

3.86194 

7.3266 

0 

-0.389316, 

-2.83213 

0 

0 

0 

i = 3 

- 1.386 

9.10716 

34.8581 

28.8994 

-5.16902 

-18.3751 

- 28.0853 

0.746189 

5.42825 


Table 8 . Numerical values of the coefficients for Eq. (2.2) for the vector current with n/ = 5. 
(Standard fixed-order expansion). 
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0 

0 

0 

i = 2 

0.999196 

1.74376 

1.19914 

0 

-0.956309 

- 1.57601 

0 

0 

0 

i = 3 

-3.19148 

1.49991 

6.92794 

1.93195 

-5.03603 

-8.67158 

-4.5967 

1.83292 

3.02069 

n = 3 

i = 0 

0.804393 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.257044 

0.804393 

0 

0 

0 

0 

0 

0 

0 

i = 2 

0.817928 

1.4748 

1.17307 

0 

- 0.492667 

- 1.54175 

0 

0 

0 

i = 3 

-1.97558 

0.0722677 

6.44038 

1.88995 

-3.75658 

- 7.59737 

-4.49678 

0.944279 

2.95503 

n = 4 

i = 0 

0.809673 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.111301 

0.809673 

0 

0 

0 

0 

0 

0 

0 

i = 2 

0.615164 

1.33706 

1.18077 

0 

-0.213327 

- 1.55187 

0 

0 

0 

i = 3 

-1.36612 

- 0.923734 

6.26766 

1.90236 

-2.62711 

- 7.08209 

- 4.5263 

0.408876 

2.97442 


Table 9. Numerical values of the coefficients for Eq. (2.3) for the vector current with Uf = 5. 
(Linearized expansion). 



[Cv]Z 

\Cv\Z 






\Cvtl 


n = 1 

i = 0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

1.19753 

1 

0 

0 

0 

0 

0 

0 

0 

i = 2 

3.1588 

4.71142 

1.45833 

0 

-2.29527 

-1.91667 

0 

0 

0 

i = 3 

1.32698 

14.7867 

13.2036 

2.34954 

-15.0028 

- 20.4771 

-5.59028 

4.39926 

3.67361 

n = 2 

i = 0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.60679 

1 

0 

0 

0 

0 

0 

0 

0 

i = 2 

2.42875 

4.12068 

1.45833 

0 

-1.16301 

-1.91667 

0 

0 

0 

i = 3 

1.7702 

11.9808 

12.3421 

2.34954 

-10.7766 

- 18.2126 

-5.59028 

2.22911 

3.67361 

n = Z 

i = 0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.31955 

1 

0 

0 

0 

0 

0 

0 

0 

i = 2 

1.65593 

3.83344 

1.45833 

0 

-0.612471 

- 1.91667 

0 

0 

0 

i = 3 

1.53408 

10.1987 

11.9232 

2.34954 

-7.11997 

- 17.1115 

- 5.59028 

1.1739 

3.67361 

n = 4 

i = 0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.137464 

1 

0 

0 

0 

0 

0 

0 

0 

i = 2 

1.0347 

3.65135 

1.45833 

0 

- 0.263473 

- 1.91667 

0 

0 

0 

i = 3 

0.744809 

8.93759 

11.6576 

2.34954 

- 4.29854 

- 16.4135 

- 5.59028 

0.504989 

3.67361 


Table 10. Numerical values of the coefficients for Eq. (2.4) for the vector current with n/ = 5. 
(Iterative linearized expansion). 
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[Canin’.- 



[tapin’.' 

[Cp]l^ 

\CptA 

\Cptl 

n = 1 

i = 0 

0.5333 

0 

0 

0 

0 

0 

0 

0 

0 

2=1 

2.0642 

1.06667 

0 

0 

0 

0 

0 

0 

0 

i = 2 

7.23618 

5.89136 

2.17778 

0 

- 4.30041 

-2.22222 

0 

0 

0 

i = 3 

7.06593 

29.1882 

19.5609 

4.47654 

- 36.7734 

- 27.9695 

- 9.07407 

8.95919 

4.62963 

n = 2 

i = 0 

0.30476 

0 

0 

0 

0 

0 

0 

0 

0 

2=1 

1.21171 

1.21905 

0 

0 

0 

0 

0 

0 

0 

i = 2 

5.9992 

6.86166 

3.70794 

0 

-2.5244 

-2.53968 

0 

0 

0 

i = 3 

14.5789 

36.2468 

31.4945 

10.0938 

-28.8842 

-32.5014 

-15.4497 

5.25916 

5.29101 

T2 = 3 

i = 0 

0.20318 

0 

0 

0 

0 

0 

0 

0 

0 

2=1 

0.71276 

1.21905 

0 

0 

0 

0 

0 

0 

0 

i = 2 

4.26701 

6.29135 

4.92698 

0 

- 1.48491 

-2.53968 

0 

0 

0 

i = 3 

13.3278 

34.8305 

38.066 

16.697 

-20.066 

-30.1251 

- 20.5291 

3.09356 

5.29101 

n = 4 

i = 0 

0.147763 

0 

0 

0 

0 

0 

0 

0 

0 

2=1 

0.401317 

1.18211 

0 

0 

0 

0 

0 

0 

0 

i = 2 

2.91493 

5.1643 

5.95979 

0 

- 0.836077 

-2.46272 

0 

0 

0 

i = 3 

9.9948 

29.5129 

40.2459 

24.1703 

-13.4331 

-25.3105 

-24.8325 

1.74183 

5.13067 


Table 11. Numerical values of the coefficients for Eq. (2.2) for the pseudoscalar correlator with 
Uf = 4. (Standard fixed-order expansion). 



[Cp]Z 

[Cp]Z 

[Cp]^’° 

\Cp]Z 

[Cp]°’,( 

[Cpln’,- 

Ppiy 

\Cptl 


n = 1 

i = 0 

0.730297 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

1.41326 

0.730297 

0 

0 

0 

0 

0 

0 

0 

i = 2 

3.58681 

2.62028 

1.12587 

0 

- 2.94429 

- 1.52145 

0 

0 

0 

i = 3 

-2.10344 

11.3262 

8.59338 

1.93901 

-19.4793 

- 13.2609 

-4.69114 

6.13394 

3.16969 

n = 2 

i = 0 

0.743002 

0 

0 

0 

0 

0 

0 

0 

0 

1 = 1 

0.738531 

0.743002 

0 

0 

0 

0 

0 

0 

0 

i = 2 

2.55535 

1.96655 

1.14546 

0 

-1.53861 

- 1.54792 

0 

0 

0 

i = 3 

0.536147 

6.35978 

7.66478 

1.97274 

-13.0167 

- 10.5778 

-4.77276 

3.20543 

3.22484 

n = 3 

i = 0 

0.766734 

0 

0 

0 

0 

0 

0 

0 

0 

2=1 

0.448296 

0.766734 

0 

0 

0 

0 

0 

0 

0 

i = 2 

2.02851 

1.71554 

1.18205 

0 

- 0.93395 

- 1.59736 

0 

0 

0 

i = 3 

1.94168 

4.12818 

7.42578 

2.03575 

-9.89039 

- 9.60801 

-4.9252 

1.94573 

3.32784 

n = 4 

i = 0 

0.787401 

0 

0 

0 

0 

0 

0 

0 

0 

2=1 

0.267317 

0.787401 

0 

0 

0 

0 

0 

0 

0 

i = 2 

1.624 

1.56872 

1.21391 

0 

- 0.55691 

- 1.64042 

0 

0 

0 

i = 3 

2.58251 

2.65675 

7.3283 

2.09062 

- 7.62431 

-9.06256 

-5.05796 

1.16023 

3.41754 


Table 12. Numerical values of the coefficients for Eq. (2.3) for the psendoscalar current with 
n/ = 4. (Linearized expansion). 
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[<5p]n’.- 

[Cp]°’.'. 


[Cp]^’] 

\Cptl 


n = 1 

i = 0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

1.93519 

1 

0 

0 

0 

0 

0 

0 

0 

i = 2 

8.78182 

5.58796 

1.54167 

0 

-4.03164 

- 2.08333 

0 

0 

0 

i = 3 

9.22126 

25.7977 

15.8503 

2.65509 

-42.7996 

- 26.4915 

- 6.42361 

8.39924 

4.34028 

n = 2 

i = 0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.99398 

1 

0 

0 

0 

0 

0 

0 

0 

i = 2 

5.42718 

4.64676 

1.54167 

0 

-2.07079 

-2.08333 

0 

0 

0 

i = 3 

11.733 

19.005 

14.3993 

2.65509 

-25.8023 

-22.5698 

-6.42361 

4.31416 

4.34028 

n = 3 

i = 0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.584683 

1 

0 

0 

0 

0 

0 

0 

0 

i = 2 

3.81501 

4.23746 

1.54167 

0 

-1.21809 

- 2.08333 

0 

0 

0 

i = 3 

11.0126 

15.8978 

13.7683 

2.65509 

-17.7717 

- 20.8644 

-6.42361 

2.53768 

4.34028 

n = 4 

i = 0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.339493 

1 

0 

0 

0 

0 

0 

0 

0 

i = 2 

2.74147 

3.99227 

1.54167 

0 

- 0.707277 

- 2.08333 

0 

0 

0 

i = 3 

9.51996 

13.9286 

13.3903 

2.65509 

- 12.512 

- 19.8428 

-6.42361 

1.47349 

4.34028 


Table 13. Numerical values of the coefficients for Eq. (2.4) for the pseudoscalar current with 
Uf = A. (Iterative linearized expansion). 





\Rv?r:° 


[Rvt^ 



[Rvt^i 


n = 1 

i = 0 

0.428571 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.0137566 

0.857143 

0 

0 

0 

0 

0 

0 

0 


1.56736 

1.44418 

1.75 

0 

- 0.0286596 

- 1.78571 

0 

0 

0 


1.72775 

1.32513 

1.67857 

0 

- 0.0263668 

- 1.64286 

0 

0 

0 


- 4.79526 

2.66831 

9.00161 

3.59722 

-6.57481 

-8.76742 

-7.29167 

0.0597075 

3.72024 

L — vj 

- 3.53855 

1.29072 

7.81479 

3.26389 

- 6.65629 

-7.1511 

-6.43452 

0.0505364 

3.14881 

n = 2 

i = 0 

0.592593 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

-0.302139 

1.18519 

0 

0 

0 

0 

0 

0 

0 


0.718416 

1.35457 

2.41975 

0 

0.629455 

- 2.46914 

0 

0 

0 


1.06685 

1.18996 

2.32099 

0 

0.579099 

-2.2716 

0 

0 

0 


- 1.59081 

- 1.60786 

11.1353 

4.97394 

- 2.02404 

- 9.44651 

- 10.0823 

- 1.31137 

5.14403 

L - ij 

0.404636 

- 3.06309 

9.54776 

4.51303 

- 3.35942 

- 7.42572 

-8.89712 

- 1.10994 

4.35391 

n = S 

i = 0 

0.681818 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

- 0.557447 

1.36364 

0 

0 

0 

0 

0 

0 

0 


-0.119176 

1.13889 

2.78409 

0 

1.16135 

- 2.84091 

0 

0 

0 


0.369655 

0.949499 

2.67045 

0 

1.06844 

- 2.61364 

0 

0 

0 


- 1.61506 

- 5.30928 

11.9551 

5.72285 

2.28504 

-9.12039 

- 11.6004 

- 2.41948 

5.91856 

i - O 

1.3858 

- 6.67953 

10.1636 

5.19255 

- 0.0698462 

-6.9352 

- 10.2367 

- 2.04785 

5.00947 


Table 14. Numerical values for the coefficients of the standard fixed-order expansion of the ratios 
of vector moments. We display results for the vector current with ny = 4, (5) for the upper (lower) 
number. 
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\Rp\Z 

[Rp]Z 

\Rpt" 

\Rptf!. 


\Rpt^ 

tapin’.- 

[Rptl 

[Rptl 

n = 1 

i = 0 

0.571429 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.0603175 

1.14286 

0 

0 

0 

0 

0 

0 

0 

i = 2 

3.262 

2.00952 

2.33333 

0 

-0.125661 

- 2.38095 

0 

0 

0 

i = 3 

6.32118 

6.21577 

12.1735 

4.7963 

- 13.7852 

- 12.0397 

-9.72222 

0.261795 

4.96032 

n = 2 

i = 0 

0.666667 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

-0.311887 

1.33333 

0 

0 

0 

0 

0 

0 

0 

i = 2 

2.1179 

1.57993 

2.72222 

0 

0.649765 

- 2.77778 

0 

0 

0 

i = 3 

9.55945 

1.01989 

12.6416 

5.59568 

- 7.82395 

- 10.8608 

-11.3426 

- 1.35368 

5.78704 

n = 3 

i = 0 

0.727273 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

-0.57611 

1.45455 

0 

0 

0 

0 

0 

0 

0 

i = 2 

1.09403 

1.25182 

2.9697 

0 

1.20023 

-3.0303 

0 

0 

0 

i = 3 

9.74702 

-3.08269 

12.8277 

6.10438 

-2.71009 

- 9.88258 

-12.3737 

- 2.50048 

6.31313 


Table 15. Numerical values for the coefficients of the linearized expansion of the ratios of pseu¬ 
doscalar moments with n/ = 4. 



[Rvt^i 

livin’.' 


livin’.' 


\Rvtl. 




n = 1 

i = 0 

0.654654 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.0105068 

0.654654 

0 

0 

0 

0 

0 

0 

0 

A _ 9 

1.19701 

1.0925 

1.00926 

0 

-0.0218891 

- 1.36386 

0 

0 

0 

Z — ^ 

1.31951 

1.00158 

0.954703 

0 

- 0.020138 

- 1.25475 

0 

0 

0 

-j — 

-3.68165 

0.823415 

5.76639 

1.73817 

-5.02124 

-6.65245 

- 4.20524 

0.0456024 

2.84138 

L — O 

- 2.72379 

- 0.349776 

4.95174 

1.53813 

-5.0835 

-5.42147 

-3.6597 

0.0385979 

2.40494 

n = 2 

i = 0 

0.7698 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

-0.196245 

0.7698 

0 

0 

0 

0 

0 

0 

0 

A _ 9 

0.441611 

1.07606 

1.18678 

0 

0.408843 

- 1.60375 

0 

0 

0 

Z — ^ 

0.667925 

0.969147 

1.12263 

0 

0.376136 

- 1.47545 

0 

0 

0 


- 0.92068 

- 1.21163 

6.45905 

2.04389 

- 1.21043 

-6.95338 

- 4.9449 

-0.851757 

3.34115 

1 — o 

0.433093 

- 2.4104 

5.5185 

1.80867 

-2.08612 

-5.57542 

- 4.3034 

-0.720927 

2.82795 

n = 3 

i = 0 

0.825723 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

-0.337551 

0.825723 

0 

0 

0 

0 

0 

0 

0 

A _ 9 

-0.141159 

1.02719 

1.27299 

0 

0.703232 

- 1.72026 

0 

0 

0 

Z — ^ 

0.154843 

0.912501 

1.20418 

0 

0.646973 

- 1.58264 

0 

0 

0 


- 1.03567 

-2.65386 

6.7324 

2.19237 

1.67114 

-6.92913 

-5.30412 

- 1.46507 

3.58387 

L — O 

0.902443 

-3.82647 

5.73411 

1.94007 

0.222185 

- 5.49342 

-4.61602 

- 1.24003 

3.03338 


Table 16. Numerical values for the coefficients of the linearized expansion of the ratios of vector 
moments. We display results for the vector current with n/ = 4, (5) for the upper (lower) number. 









\Rptf:. 


\Rptl 

\Rpt^. 

[Rpt^. 

[Rpt^^ 

n = 1 

i = 0 

0.7559290 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.0398962 

0.755929 

0 

0 

0 

0 

0 

0 

0 

i = 2 

2.15656 

1.28928 

1.16539 

0 

-0.0831172 

- 1.57485 

0 

0 

0 

i = 3 

4.06725 

1.88674 

6.70126 

2.00706 

-9.11366 

- 7.79727 

- 4.85579 

0.173161 

3.28094 

n = 2 

i = 0 

0.816497 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

-0.190991 

0.816497 

0 

0 

0 

0 

0 

0 

0 

i = 2 

1.27461 

1.1585 

1.25877 

0 

0.397898 

- 1.70103 

0 

0 

0 

i = 3 

6.15209 

- 0.379067 

6.87731 

2.16787 

-4.6981 

- 7.44666 

- 5.24486 

- 0.828954 

3.54382 

n = 3 

i = 0 

0.852803 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

- 0.337775 

0.852803 

0 

0 

0 

0 

0 

0 

0 

i = 2 

0.574538 

1.07172 

1.31474 

0 

0.703697 

- 1.77667 

0 

0 

0 

i = 3 

5.94225 

- 1.95744 

6.96991 

2.26427 

- 1.31021 

- 7.20157 

-5.47807 

- 1.46604 

3.7014 


Table 17. Numerical values for the coefficients of the linearized expansion of the ratios of pseu¬ 
doscalar moments with n/ = 4. 



\Rvt^i 


[Rvt^i 


\Rvt^i 

\Rvt\ 

\Rvtl 


\RvtA 

n = 1 

i = 0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.0160494 

1 

0 

0 

0 

0 

0 

0 

0 

i = 2 

1.86056 

3.66883 

1.54167 

0 

- 0.0334362 

- 2.08333 

0 

0 

0 


2.04768 

3.52994 

1.45833 

0 

- 0.0307613 

- 1.91667 

0 

0 

0 


- 1.85046 

11.8662 

12.8916 

2.65509 

-7.80382 

-18.4951 

-6.42361 

0.0696588 

4.34028 

L — O 

- 0.0174324 

9.52394 

11.4806 

2.34954 

- 7.88822 

- 15.9481 

-5.59028 

0.0589592 

3.67361 

n = 2 

i = 0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

- 0.254929 

1 

0 

0 

0 

0 

0 

0 

0 

i = 2 

0.0638103 

3.39785 

1.54167 

0 

0.531103 

-2.08333 

0 

0 

0 


0.357801 

3.25896 

1.45833 

0 

0.488615 

- 1.91667 

0 

0 

0 

7 — 

-2.11686 

9.07965 

12.4739 

2.65509 

0.552022 

-17.366 

-6.42361 

- 1.10646 

4.34028 

1 — o 

0.322201 

6.88187 

11.0854 

2.34954 

- 0.755491 

- 14.9093 

-5.59028 

-0.936512 

3.67361 

n = 3 

i = 0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

- 0.408795 

1 

0 

0 

0 

0 

0 

0 

0 

i = 2 

- 0.988542 

3.24398 

1.54167 

0 

0.851656 

- 2.08333 

0 

0 

0 


- 0.630066 

3.10509 

1.45833 

0 

0.783523 

- 1.91667 

0 

0 

0 

7 — 

-5.11183 

7.46526 

12.2367 

2.65509 

5.43048 

- 16.7249 

- 6.42361 

- 1.77428 

4.34028 

L — O 

- 1.87845 

5.35334 

10.861 

2.34954 

3.40317 

- 14.3195 

- 5.59028 

- 1.50175 

3.67361 


Table 18. Numerical values for the coefficients of the iterative linearized expansion of the ratios 
of vector moments. We display results for the vector current with ny = 4, (5) for the upper (lower) 
number. 
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[Rpfn!. 


IRpfrl!. 

[«p1°’,- 

[RpU 

n = 1 

i = 0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

0.0527778 

1 

0 

0 

0 

0 

0 

0 

0 

i = 2 

2.95841 

3.70556 

1.54167 

0 

-0.109954 

-2.08333 

0 

0 

0 

i = 3 

11.4629 

13.0982 

12.9483 

2.65509 

- 12.4961 

- 18.6481 

-6.42361 

0.22907 

4.34028 

n = 2 

i = 0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

-0.233915 

1 

0 

0 

0 

0 

0 

0 

0 

i = 2 

1.09324 

3.41886 

1.54167 

0 

0.487324 

-2.08333 

0 

0 

0 

i = 3 

8.77473 

10.1858 

12.5063 

2.65509 

- 3.80468 

- 17.4536 

-6.42361 

- 1.01526 

4.34028 

n = 3 

i = 0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

i = 1 

- 0.3960760 

1 

0 

0 

0 

0 

0 

0 

0 

i = 2 

-0.118446 

3.2567 

1.54167 

0 

0.825158 

- 2.08333 

0 

0 

0 

i = 3 

4.92499 

8.38182 

12.2563 

2.65509 

1.76427 

- 16.7779 

-6.42361 

- 1.71908 

4.34028 


Table 19. Numerical values for the coefficients of the iterative linearized expansion of the ratios 
of pseudoscalar moments with n/ = 4. 
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